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Abstract 


This report presents the theory, results and user instruc- 
tions for the aerodynamic program. The theory is based on linear 
^fting surface theory and the method is the kernel function 
The program is applicable to multiple interfering surfaces which 
may be coplanar or non-cop lanar. Local linearization is used to 

problems without shocks. For cases with 
imbedded shocks, the appropriate boundary conditions are added 
to account for the flow discontinuities. The data describing 
non-uniform flow fields must be input from some other source 
such as experiment or a finite difference solution. The results 
n j small linear perturbations about non-linear 

flow fields. The method is applied to a wide variety of problems 
for which It IS demonstrated. to be significantly superior to the 
uniform flow method. The program user instructions are given in 
the last appendix for easy access. 
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A STEADY AND OSCILLATORY KERNEL FUNCTION 


METHOD FOR INTERFERING SURFACES IN SUBSONIC, 
TRANSONIC AND SUPERSONIC FLOW 


By Atlee M. Cxinningham, Jr. 

Fort Worth Division of General Dynamics 


SUMMARY 


This report presents the theory, results, and user instruc- 
tions for a computer program to calculate steady and unsteady 
aerod3mamics on interfering surfaces in subsonic, transonic or 
supersonic flow. The theory is based on linear lifting surface 
theory and the method is the kernel function. The program is 
applicable to multiple interfering surfaces which may be coplanar 
or non-cop lanar. Local linearization is used to treat non-uniform 
flow problems without shocks. For cases with imbedded shocks, 
the appropriate boundary conditions are added to account for the 
flow discontinuities. The data describing non-uniform flow fields 
must be input from some other source such as experiment or a 
finite-difference solution. The results are in the form of small 
linear perturbations about non-linear flow fields. The method is 
applied to a wide variety of problems for which it is consistently 
demonstrated to be significantly superior to the uniform flow 
method. 

All equations used in the method are sximmarized in the 
appendices along with detailed derivations. A set of instruc- 
tions for using the program is presented in the last appendix 
for easy access. The instructions are included in this report 
so that reference to pressure function or surface types and 
example cases can be easily made. 

In addition to the non-uniform flow capability, the program 
can be used to obtain uniform flow solutions with a single option 
change. Interference effects can be calculated for coplanar, 
non-coplanar and intersecting planar surfaces. The aerod3mamic 
input geometry data format permits the input of arbitrarily 
arrayed surfaces which are constrained only to being streamwise. 
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Each surface may or may not have an image surface, hence, 
as)mHnetric configurations may be constructed. 

The program has a built-in interpolation scheme for struc- 
tural mode shapes. The scheme uses a surface spline fit over 
various structural surfaces or regions. This will permit inter- 
polation of modes for total airplane configurations with all- 
movable and vertical surfaces. 

Generalized forces can be calculated in unsteady flow and 
pressure distributions can be obtained in both steady and un- 
steady flow. The solutions obtained are essentially independent 
in computer cost of the number of modes or downwash vectors input. 
Once the aerod}mamic matrices are conq>uted, inverted and saved on 
a magnetic tape, they can be used on subsequent problems for very 
little cost as long as Mach niimber, reduced frequencies, and aero- 
dynamic geometry remain unchanged. Thus, the method is tailored 
for design applications where the structural mode shapes change 
continually for structural changes and payload variations while 
aerod 5 mamic parameters remain constant. 


INTRODUCTION 


In recent years, interest has grown considerably in the 
desire to fly efficiently in the high subsonic regime. As a 
result, the need has increased for better unsteady transonic 
aerod 5 mamic tools so that flutter and djmamic response charac- 
teristics can be more accurately predicted in this flow regime. 
Presently, these characteristics are predicted only with methods 
which are based on linearized theory in uniform potential flow. 
In addition, since buffet and limit cycle flutter appear to be 
similar in experimental flutter and buffet testing, it is impor- 
tant that their distinction be better understood. 

The characteristic of transonic flow which causes the 
greatest difficulty when attempting to apply tmiform flow theory 
to such problems is the presence of shocks imbedded in the flow. 
Such a gradient in velocity as that which exists across a shock 
is no longer small, thus, linear theory methods (ref. 1) cannot 
account for this phenomenon and hence become invalid. Finite 
difference methods or other iterative schemes can accoxmt for 
such discontinuities but they are usually very expensive to use 
in terms of computer time required (refs. 2, 3, 4 and 5). More- 
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over, if they are used for flutter or dynamic response analyses 
where solutions must be computed for 10 or more frequencies, the 
computer costs quickly become astronomical. 

For example, a simple cantilevered wing flutter analysis 
with four natural modes and 10 frequencies could require 20 or 
more hours of computer time on an IBM 360/65 for a single flutter 
solution. Thus, such an approach is not well suited for solving 
unsteady transonic aerodynamic problems in a practical sense, 

A study was conducted to investigate the feasibility of 
using combined subsonic and supersonic linear theory as a means 
for solving unsteady transonic flow problems economically 
(ref, 6), In the method developed, a wing over which the flow 
was mixed supersonic and subsonic with imbedded shocks was 
treated as an array of general aerod)niamic lifting surface 
(GALS) elements. Each element was allowed to have mutual inter- 
ference with the other elements. Also, each was assigned a 
different Mach number, either subsonic or supersonic, and its 
downwash was modified accordingly. The Mach number distribution 
and shock geometry was obtained from either experiment or a 
finite difference solution, hence the method was used to predict 
unsteady perturbations about known steady mean flows. Once 
assembled, the solution proceeded in a manner identical to ordi- 
nary aerod 3 mamic interference methods (ref, 7), The frequency 
sweep could be performed at about the usual cost of a standard 
subsonic or supersonic unsteady aerod 5 mamic analysis which is 
less than one hour - usually about 10 minutes - as opposed to 
20 hours or more for a finite difference solution. As a result 
of the feasibility study the computer procedure documented in 
this report was developed. 

The theory, results and program utilization are given in 
this volxime for the subsonic, transonic and supersonic aero- 
d 3 mamic program. The program is described in Volxjme II of this 
report. This program is applicable to steady and tmsteady flow 
over multiple arbitrarily arrayed lifting surfaces. The Mach 
range is subsonic, mixed transonic with or without shocks, and 
supersonic. In addition to uniform flow subsonic and supersonic 
solutions, non-uniform Mach number distributions can be used 
which lead to significant improvements in solution accuracy. 

The method is based on a kernel function technique which uses 
asstmied pressure functions with unknown coefficients. By match- 
ing flow tangency boundary conditions at the control points, the 
unknown coefficients are uniquely determined. The Mach number 
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and reduced frequency at the control point are used to calculate 
the aerodynamic influence coefficients at that point. The 
presence of a normal shock is simulated by a line doublet which 
represents the load induced by shock movement. The appropriate 
steady or unsteady normal shock boundary conditions, as derived 
in this report, are satisfied across the shock along the surface 
of the wing. In the application of the transonic method to 
several cases, the solutions are shown to be significantly 
superior to tmiform flow theory solutions with a relatively 
small increase in computer costs. 


SYMBOLS 


free stream speed of soimd, meters /second 
conical coordinate 


AR 


aspect ratio 


^REF 

b(^) 


reference length, meters - usually % wing chord 
for 2-dimensional flow or % MAC for finite wings 
in 3 -dimensional flow 

wing semi-chord at span station 5f, 
non-dimensionalized . by bj^ 


pressure coefficient 


b(S,y) 


i = 



mode amplitude at point x,y 
non-dimensionalized by b^Ep 


k = 


6>bREF 

U 


m 


tan 


reduced frequency 


M 


= U 


a 


Mach number 


4 


MAC 


mean aerodynamic chord, meters 


PqCr,^) 

p 

A Pq( £, 2? 
q = 

2 

r = J(y-’)7 

So 

U 

u,v,w 

Wp(sr,3r) 

Wp(Sc,5r) = 
x,y,z 

/w 

x,y,z 

^o »yo »^o 

A Cp 
a 

iS =/T^ 
/3^ = 1-m2 


supersonic pressure weighting function 

in the plane of the qth surface, non-dimensional 

o 

pressure. Newtons /meter‘d 

) lifting pressure amplitude in the plane 
of the qth surface. Newtons /meter^ 

dynamic pressure. Newtons /me ter^ 


)2 


+ (z- C 



h 


non- dimens ionalized by 


wing semi-span, non-dimensionalized by bRgF 


free stream velocity, meters/second 

velocity components in the x,y,z directions, 
respectively, meters /second 

amplitude of the oscillatory downwash 
normal to the pth surface, meters /second 

Wp(x,y) /U 

cartesian coordinate location of the downwash 
point in the kernel function (x is in the 
direction of U) , non-dimensionalized by b^Ep 

coordinates in the plane of the p^^ surface 
with z perpendicular (see fig. 3.) 

distance from an influence point to the downwash 
point, (x- ^ ) , (y- 17 ) , (z- O » non-dimensional 

lifting pressure coefficient 


angle of attack, degrees 
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Ale» A.j,g leading and trailing edge sweep angles, degrees 


^ , 1? , C location of an influence (or integration) point 

in the kernel function, non-dimens ionalized by 

i, X coordinates in the plane of the q^h surface 

with X perpendicular (see fig. 3.) 




location of the mid-chord at span station 
non-dimensionalized by bj^p 


i. 2?. i 

8 

U) 


( - (.(v) ± X 

b(v) ’ So’ So 

limits of the spanwise integration as determined 
by geometry and the Mach hyperbola 

chordwise variable of integration, non-dimensional 

rotational frequency, radians /second 


Subscripts : 


L 

LE 

MC 

P 

q 

TE 

xy 


local value 
leading edge 

Mach cone (or hyperbola) boundary 

downwash surface 

integration surface 

trailing edge 

local value at point x,y 

free stream conditions 
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TRANSONIC FLUTTER AND DYNAMIC RESPONSE ANALYSIS 


Before proceeding with describing the transonic method, 
it is appropriate that a brief discussion be given concerning 
the problems associated with transonic flutter and dynamic 
response analysis. 

In conventional d 5 mamic analysis, the points to be considered 
are defined over a Mach number-altitude envelope. Since the 
analysis methods use linear theory aerodynamics, there is no 
coupling between the unsteady and steady flow fields, and altitude 
is accounted for only through the variation of air density. In 
transonic flow this is no longer true. 

The prediction of unsteady pressure distributions induced by 
a surface oscillating in a mixed transonic flow is complicated by 
the strong coupling between the steady and unsteady flow fields. 
The steady flow fields are in turn drastically modified by Mach 
number, altitude, thickness, camber, twist, angle of attack, 
planform geometry, interference, and static aeroelastic and 
boundary layer effects. For a given Mach-altitude point and a 
fixed configuration, however, the only additional variable is the 
angle of attack, a. Thus, transonic flutter and d 3 mamic response 
analysis must be performed over a three-dimensional envelope as 
specified by Mach-altitude- a conditions. 

Transonic analysis is further complicated by the need to 
compute unsteady pressures on the entire configuration, that is 
both upper and lower surfaces, in most cases. Through changes in 
the static aeroelastic deformation and boundary layer, altitude 
effects in the shock structure and flow fields are as significant 
as those due to Mach number and a . Thus , a new mean flow field 
is needed for each Mach-altitude- a condition. 

The use of finite difference or other iterative schemes such 
as those given in references 2, 3, 4 or 5 would pose an obvious 
solution to this problem. Assuming that an unsteady verision 
could be developed for finite wings, it might be assumed that 
about 30 to 60 minutes on an IBM 360/65 would be required for a 
single aerod 5 niamic solution (ref. 4). A single solution would, 
however, refer to one Mach-altitude- a -frequency-mode condition. 
Thus for a single 10 frequency, 4 mode flutter solution for a 
simple cantilever wing, 20 to 40 hours of computer time would be 
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required. Expanding this to two Mach n^Jmbers, two altitudes, and 
two a's, the total cost would be 160 to 320 hours of computer 
time. This cost is on the order of magnitude of a flutter model 
test program and the same cost would be encountered for redesign 
evaluations. Since one of the primary objectives for using theo- 
retical flutter methods is to reduce costs by minimizing the 
requirement of flutter model tests, then the use of the finite 
difference approach would defeat its intended purpose. 

Although linear theory methods cannot solve the highly non- 
linear mean flow problem, they can solve the small perturbation 
problem about the mean flow. Thus, a more realistic and yet 
economical approach would be to use the elaborate schemes to 
predict the mean flow fields and linear "transonic” theory to 
perform the flutter analysis at each Mach-altitude- a point. 

The cost of using linear theory methods is primarily a function 
of the number of frequencies and is practically independent of 
the number of modes. Also, once the aerodynamic matrices are 
calculated, they may be used repeatedly to evaliiate design changes 
at very little cost. 

Returning to the cantilevered wing, for a linear theory 
solution of about 1 minute per frequency, the cost for flutter 
analysis would be 10 minutes per Mach-altitude- a point. For the 
total, 4 to 8 hours would be required for the steady mean flow 
calculations and 1.33 hours for flutter analyses. Comparing 5.33 
to 9.33 hours against 160 to 320. hours, it is clear that the 
hybrid approach would be economically very attractive. The re- 
sults in this report will also demonstrate that the use of linear 
theory would not degrade the accuracy. 


THE TRANSONIC KERNEL FUNCTION METHOD 


The fundamental problem to be treated in this report is the 
development of a technique to solve the integral equation that 
relates the normal velocity imposed by boundary conditions with 
the load distribution on an arbitrary array of planar lifting 
surfaces in a subsonic, supersonic or mixed transonic flow. 

The equation may be written as 


Wp(x,y,z) 

U 


1 Q 



C )K(x-^,y-/7,2-^,k,M)ds 


( 1 ) 
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for Q total surfaces. The downwash Wp(x,y,z) is the velocity 
normal to the lifting surface at control point (x,y,z). 

The function Apq( ^ yV , is the normal lift distribution on 
the qth lifting surface at load (or integration) point ( ^ ^ , C ) • 
The kernel function K( ) is the influence function which is 
actually the velocity field due to an elemental normal load at 
p^int ( ^ , V , C ) £n the q^h surface. The unknown quantity is 
^P( f » ■’7 > O 9-iid Wp(x,y,z) is prescribed by the boundary condi- 
tions on the lifting surface due to surface slope and motion. 

The method used to solve equation (1) in this report is based 
on a collocation technique. The unknown pressure function is 
assumed to be composed of a series of pol)momials weighted by a 
user selected weighting function that is characteristic of each 

surface. The non-planar kernel function is used for which 
the Mach number and reduced frequency are determined by those 
values at the downwash control point, x,y,z. The necessary equa- 
tions for evaluating equation (1) are given in the Appendix. 

The algorithm for linking subsonic and supersonic linear 
theory solutions together is based on the following two assump- 
tions : 

1. The appropriate Itoch number for computing downwash 
at a point is the Mach number of that point. 

2. The reduced frequency, k^^ = ^ is modified 

according to the local velocity such that co 
is held constant. This is approximated as 



The first assumption is justified by the fact that, for any 
given pressure distribution, the integrated kernel function- 
pressure function product rapidly becomes independent of Mach 
number as distance increases either upstream or downstream from 
a loaded region. The second assumption is mandatory since the 
physical frequency, o», must be held constant. The use of the 
Mach number ratio rather than velocity ratio will result in a 
small overestimate (typically about 3% to 4%) of the effect of 
velocity change. 
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With the two basic assumptions, the computational algorithm 
becomes a simple problem of testing the Mach number of the down- 
wash point. If the downwash point is supersonic, then the self- 
induced downwash as well as all interference effects at that point 
are computed with the supersonic kernel function regardless of the 
interfering surface's Mach number. Likewise, if the downwash 
point is subsonic, the subsonic kernel function at that Mach 
number is used. The value of k in the kernel function is also 
determined by the local Mach number since the downwash surface 
sees the same value of w regardless of what surface the distur- 
bance is emitted from. 


The correct form of the downwash-p res sure function integral 
equation which embodies the above assumptions can be derived from 
the nonlinear partial differential equation that governs the flow 
potential. The equation in vector form is 


1 d^<t> 



where 


a 


2 





[(q • V)q 


y-1 

~r 



ao = stagnation speed of sound 


q = V ^ 

V = - ± + - ± + t ^ 

-^y dy dz 

Let 


<P 


ti>, + (i) , 




1 


and 


Uvv 

M, = - -Jii * M 

1 a xy 

i xy 

lx xy 


Defined at the 
downwash point 


where Uxy is the mean flow field velocity parallel to the lifting 
surface at the downwash point location, (x,y). The first of the 
two basic assimiptions imples that the derivatives of the mean 
flow field velocity components are small relative to velocity 
component i.e. 
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( 2 ) 



■^Ixy’ '*’1x2^ ~ (‘'■lyx' *lyy' ‘*'lyz> 
'^'^Izx- V’ *lzz>~^K« <^1^ 


Employing equation ( 2 ), it is possible to reduce the nonlinear 
equation to the following first order nondimensional form for (f> 


V^<p - <t> 


- nkiKi <f>^ = 0 


where 


<#>» 


0 


^Ix bREF 


V- bREp-V 


<^x = 




It should be noted that all equations have been defined at the 
downwash point. 


Next, the acceleration potential equations are developed 
also at the downwash point with one exception. The acceleration 
potential for harmonic motion is defined as 


(f) 

ix X 

It is desired to non-dimens ionalize 'I' with a flow field variable 
at the integration point since it is directly related to the 
pressure difference. Hence, 



where is the mean flow velocity component at the integration 
point, TOe governing differential equation for can thus be 
written as 
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for which the well known solution is 


where 


, A 5 1 


Ml 

Ml (=-« - H I 

^1 


4 


= 1 - Ml 


r'2 = 


(X-^) + [ (y-r,)2 + (2-0^1 


The strength factor A is determined from the relationship 




2 


AC 


EL 


Since flow can be assumed as nearly two dimensional in the 
vicinity of the downwash point, the relation 




, ^loc 

'/'= -s-T — ACp^ 


2</>. 




is approximately valid where and are the d 3 mamic pressures 

at infinity and the integration point, and ACpoo is the integration 
point lift amplitude divided by q^^. As a result, the solution to 
equation (3) is 


1 ^Ico -iki(x-^) 




^ Stt (f>ix ^ 


) /• " 

J aix 


• -oc 


dz 


^ EXP 
^1 


^*'1 ^2 


e 


-%X 


dX 


where ACp^(^,'»j, {) is the lifting coefficient at the integration 
point. Defining the downwash as 


1 d<f> d(f> 


(dz 


•hx dz dz ' 


The downwash-pressure integral equation becomes 




• K (x-^, y-Tj, z-^, ki. Ml) d^di 7 


where the kernel function is defined as 


K (x-f, y-rj, z-^, ki, Mi) 

az2 *^oe 


EXP 


iki ^ I Mix -Rj^ 


-ikiX 
e ^ dX 


The only difference between the above formulation of the kernel 
function and the classical form is that the ki and Mi vary with 
downwash point location. They are constant over the region of 
integration, however, for any given downwash point. These 
definitions are consistent with the two basic assiomptions. 

The presence of the term (<^l^/^lx) the integral in 
eqxiation (5) violates the first assumption that the downwash 
point Mach number dominates the influence functions. Likewise 
as the integration points approached the downwash point, this 
term would approach unity, therefore its influence could be 
negligible. For the purposes of the study stimmarized in this 
report, it has been set equal to unity. The effect could be 
significant, however, for integration points in a supersonic 
region nearly upstream of a subsonic downwash point. The use 
of should be investigated in future stixdies to 

determine its importance. Also, rather than use (<Alx/<Al«;) 
equation (5), . . m 

<^lx/</>loo « ^y/M 

is used which will introduce a small error as previously discussed. 


13 



As a result of the above discussion, the integral eqviation 
that is actually solved by the transonic algorithm is 


-(ft ^ • 

K(x-^, y-^, Z"f, kjjy, Mxy)d|*^^ 

( 6 ) 


where ACp ( i t V ) is the pressure difference at ( f , '>7 ) divided 
by q«o . iince the algorithm permits Mxy and kxy to be used in 
computing the kernel function, thereby permitting the kernel 
function to change at each downwash point, it is possible to 
account for leading edge regions where the flow is continuously 
accelerated from subsonic to supersonic. 

Shown in figure 1 is a flow diagram of the current transonic 
flow algorithm for unsteady flow. The same logic is applicable 
to steady flow. The key ingredient is- the common set of surface 
types and pressure function types for both subsonic and supersonic 
flow as shown in figure 2. The integration schemes are also 
compatible with two exceptions. The difference in the chordwise 
integration is that in supersonic flow, the limits are from Mach 
cone to leading edge whereas in subsonic flow they are from trail- 
ing edge to leading edge. In the spanwise integration, both 
techniques are identical. 

The pressure functions are given as a set of chordwise and 
spanwise varying pol)niomials with a weighting function which is 
constructed according to the pressure , type shown in figure 2, 

The form is 


4Pq(i,!?) Soq h(i)i(5)[goq(5)fo(i) 

+ ...J 

( 7 ) 

where all items are defined in the "Symbols” section except for 
gmq(22)> fn(£)» h(£) and i(V). The coordinates, £, v, are 
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IS-0 


IS-IS+1 


Mic <1.0 


£ 


IC-0 


I 


IC-IC+1 


I 


k'-4c« 


i 

m 


Mic 


Mjq ^1.0 


L 


COMPUTE SUBSONIC 


COMPUTE SUPERSONIC 

UNSTEADY AERO 


UNSTEADY AERO 

MATRIX BAND FOR 


MATRIX BAND FOR 

DOWNWASH POINT "IC" 


DOWNWASH POINT "IC" 

ON SURFACE "IS" FOR 


ON SURFACE "IS" FOR 

ALL SURFACES 


ALL SURFACES 


1CP( IS) “NUMBER OF 

DOWNWASH POINTS 
ON SURFACE "IS" 

NSURF-TOTAL NUMBER 
OF SURFACES 



IC< ICP(IS) 


IC-IC(IS) 

IS < NSURF 


IS=NSURF 


. COMPUTE UNSTEADY AIRLOADS, 
PRESSURE DISTRIBUTIONS OR 
GENERALIZED FORCES 


Figure 1, - Flow Diagram For Computing Unsteady Aerodynamics 
in Mixed Transonic Flow 
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TYPE 1 


TYPE 2 


TYPE 3 


SURFACE TYPES 



I 


/(5) 






4 + 3 





i 

S PANWISE LOADINGS 


I /nj 



CHORDWISE LOADINGS 


Figure 2. - Surface and Loading Types for the General 
Aerodynamic Lifting Surface Element 



defined in the transformed plane of the surfaces shown in figure 3, 
The chordwise functions, h(£), are defined as 



, Type 1 

h(i) - /w 

, Type 2 

h(i) 

\l 1 + i 

, Type 3 

h(^) - 1 

, Type 4 

in relation to figure 2. 
defined as 

The spanwise functions, 

i(2) = /T- 2^ 

, Type 1, surface 1 only 

li 

, Type 2, surfaces 1, 2, 

i(2) = - 2 

, Type 3, surfaces 2, 3 

i(s) ’ h +1 

, Type 4, surfaces 2, 3 


The chordwise functions, fn( 1 ) , are simply 
fo(i) « Uo(i) = 1 
^i(i) - Ui(i) + - 2« + 1 

^n(i) - Un(i) + Vld) 



There are mq total chordwise functions corresponding to m total 
downwash points in the chordwise direction for the qth surface. 


17 






Figure 3. - Basic Geometry and Coordinate Systems 
for Interfering Surfaces 


The spanwise functions, gmq(2i), contain the unknown coefficients 
in the form 

8mq(2) “ [^oq^o(2) + + ••• j 

where Un(2?) are the Tschebychev polynomials of the second kind, 

Uo(2) = 1 
Ui(2) = 22 


Un(2) = Un-l(2) * 2 3U„_2(2) 

There are nq total spanwise functions corresponding to nq total 
downwash chords. 

These pressure functions are used regardless of Mach number 
for subsonic or transonic solutions. For supersonic flow, these 
functions may be used or the supersonic weighting function may be 
used instead. 

If the supersonic weighting function is used, the form is 


4Pq(^2) = 4pu2p(r,2)[goq(2) fo (D 
+ 8lg(2) fi (i) + • • • j 


( 12 ) 


where P( 1^,^ ) is the weighting function. The gnq(2) and fn( £ ) 
are identical to those given for the regular function in 
equation (7). 


The supersonic weighting function is based on conical flow 
theory solutions to the lift distributions on flat swept wings 
(ref. 8 ). Some liberty has been taken to simplify the expres- 
sions and yet maintain the basic characteristics. The function 
has been developed only for simple trapezoidal wings in this 
program. The following derivations provide the equations 
necessary for wings with or without clipped tips and with sub- 
sonic or supersonic leading and trailing edges. No secondary 
reflections of the Mach lines are accounted for since they tend 
to be of second order effect and can be adequately accounted for 
in the collocation solution. 
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The basic equation for leading edge and root characteristics 
is the delta wing distribution. For a subsonic leading edge, 


m = 


tan A 


LE 


< 1 


we have 


p(r,^) = 




(13) 


where 


thus 


_ Mn-n) 

^Xi 


a ^ (^-yi) tan Ale 
™ f-Xi 


= Location of leading 
edge vertex 


For a supersonic leading edge, m >1.0, 


P(r,»7) = 1 - u(a)L(m)7l-a 




(14) 


where u(a) = 1, a<l 

= 0, a^l 


L(m) 


= 1 - 


4m 


7.0 


f 7.0 
Ll.75+l/m_ 


The term y 1,75 which appears in the expression for L(m) 

is the approximation used for the exact function 

. 4maq ^ 7otq 

PrOOT ■ /3E' (m) ^ (1.75+1/m) 


where E'(m) is the complete elliptic integral of the second kind 
of modulus (l-m^)?. 
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A tip correction is included for clipped tips. For a sub- 
sonic leading edge, m<l,0, the lift distribution behind the tip 
Mach line is constant in the streamwise direction. The amplitude, 
a function of span only, is given as 

, 1 f Fg+a- H-H-a’ ) 1 ^ 

V m ^ 

where 


(15) 


a' 


Srip-(yTip-^)^ 


= value of "a” along the tip 
Mach line at span station 


For a supersonic leading edge, m >1.0, the lift distribution 
behind the tip ^fa,ch line is given as 


P(?»^) = P(r,^) sin ^ J a' ' 

DELTA '' J 

where 

= Delta wing distribution given 
D LTA equation (14) . 


(16) 


•II 



= value of *'a” relative 
to forward wing tip. 


The functions given above in equations (15) and (16) are exact 
shapes as required by conical flow theory. 

A final correction to the delta wing distribution is the 
subsonic trailing edge term. This term is approximated as a 
multiplicative function applied to the delta plus tip term. 

The function is 


p(r,^) 


p(r,?D 

DELTA+TIP 




(17) 
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for 0.0 < a” 


^ 1,0 where 


a’” 

Xte 

^TEV 


f"%E 


(Xtev+^ -x^Ev 

2i(tan Aijg) + XjgY = x position 
X position of the trailing edge 


of trailing edge 
vertex 


thus 


a’” 

m’ 


^-(f-x 
'^(1-m') 


TEV'' £ 




tan A 


< 1 


TE 


The form of this approximation is not exactly correct, however, 
it seems to be close enough for practical purposes as experience 
has shown. 


Three examples are shown in figures 4, 5 and 6 for the 
supersonic weighting function. The magnitude of the weighting 
function is adjusted uniformly in each case so that the shapes 
can be compared with other theories or experiment. The first 
example in figure 4 is for a rectangular wing, AR=2.0, in steady 
flow at M=l,2 and a: =1,0 rad. The solid line is the weighting 
function evaluated at the span stations 22=0.1, 0.5 and 0.9 with 
the equations given in this section. The Mach lines are shown 
for clarity. The symbols are values computed by the AFFDL Mach 
box program for wing- tail configurations (ref. 9). The second 
example in figure 5 is a swept tapered wing of the standard AGARD 
wing-tail configuration (the true planform is shown). The con- 
ditions are steady flow, M“1.2 and a=1.0 rad. The solid line 
is again the weighting function. The S5nnbols are results from 
the Woodward finite element method (ref. 10). The disagreement at 
the Mach line discontinuities is due to the inability of the 

element representation to conform to such characteristics 
with a reasonable number of elements. The third example shown in 
figure 6 is a trapezoidal wing with a supersonic leading edge. 

In this case, the leading edge vertex Mach cone intersects the 
tip cone. Comparison with the Woodward method and experiment 
clearly illustrates the validity of the function. 
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Moo = 1.2 
a = 1.0 RAD 
AR = 2.0 


o AFFDL MACH BOX METHOD (REF. 9) 
SUPERSONIC WEIGHTING FUNCTION 



Figure 4, - Supersonic Weighting Function for 
a Rectangular Wing 
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o WOODWARD METHOD (REF. 10) 
_ SUPERSONIC WEIGHTING 
FUNCTION 



Figure 5. - Supersonic Weighting Function 
for a- Trapezoidal Wing with 
Subsonic Leading and Trailing Edges 
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= 2.01 o EXPERIMENT (REF. 11) 

a =6.0° SUPERSONIC WEIGHTING 

AR = 1.34 FUNGTION 

WOODWARD (REF. 10) 



Figure 6. - Supersonic Weighting Function for 
a Trapezoidal Wing with 
Supersonic Leading and Trailing Edges 
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METHODOLOGY FOR INCLUSION OF IMBEDDED SHOCKS 


This section summarizes how the presence of normal shocks in 
the flow is accounted for in the linearized model. The mean flow 
field may be non-linear but the perturbations must be small enough 
such that their non-linearities are second order effects. 


Normal Shock Boundary Conditions 
and the Shock Model 

In the calculation of unsteady aerod5mamic pressures induced 
by surface motion in mixed transonic flows, the influence of 
imbedded shocks must be accounted for by satisfying the proper 
boundary conditions across the shocks on the flow perturbation 
potential. The boundary conditions are the Rankine-Hugoniot 
relation and the equality of total potential in front of and 
behind the shock. The constraint is that the shock movement 
necessary to satisfy Rankine-Hugoniot is equal to that necessary 
to maintain equality of potential. The resulting boundary con- 
dition in steady flow as derived in Appendix D is 

K<^ (18) 


<t> - Perturbation potential 
<t>-^ - Mean flow potential with shock 


where 


“f" “f* * 

<t> - K0 = ~ u4> 

X X 


K = 


<x - 


^ = 


y-1 


y+i M2(y+i) 



The superscripts, + and refer to conditions just dovmstream 
or upstream of the shock, respectively, as shown in figure 7. 
Since the right hand side of equation (18) is determined only by 
upstream perturbations and mean flow conditions, the equation 
serves simply as an additional constraint on the downstream 
perturbations . 


Considering equation (18) for a case of downstream excitation 
only, such as a control surface deflection, it can be shown that 
the existing lifting surface representation must be changed. For 
such a case, = <^x = 0» and since the downstream region is sub- 
sonic, = 0. Thus, in order to satisfy eqxaation (18), must 
be zero; however, this is not the case as shown by experiment. 

In order for ^ 0, then 0, thus another potential must be 

added to the multiple surface model for mixed transonic flow with 
shocks • 


The finite potential </>"*■ can be achieved by placing a lifting 
line (or line doublet) between the upstream and downstream sur- 
faces along the shock location. In reality, this lifting line 
represents the lift due to shock movement which is otherwise not 
accounted for in the multiple surface model. Without the “shock 
doublet", the solution yields a peak in pressures at the leading 
edge of the subsonic surface by satisfying the flow tangency 
conditions. This peak, however, violates the shock boundary 
conditions since = 0, or else it requires that K = 0 in 
equation (18), If = 0, the shock movement, \, (shown in 
Appendix D, equation (Dll)) becomes 


X 


0 + -<f 

^Ix ^Ix 


= 0 


since = 0, Thus, for the case of downstream excitation, 
the shock does not move which is contrary to experimental 
observations. Also, for K*0, it is implied that 





0 


which is inadmissible since the shock movement would be infinite 
as given by 


X » - 


_ - ^x'*’^*^X 


0T +f4~ 
Ixx Ixx 


00 
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Figure 7. - Mean and Perturbation Flow Potentials 
About an Idealized Normal Shock 
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(shown in Appendix D, equation (D13)). The final conclusion is 
that if the constraint 

k/ = 0 

is imposed, satisfaction of the shock boundary conditions becomes 
an indeterminant problem. Thus, for the case of mixed flow with 
imbedded shocks, these boundary conditions must be satisfied in 
order to properly represent the flow fields. They can be satis- 
fied only by inclusion of the shock doublet. 

The normal shock boundary conditions for oscillatory flow 
perturbations are similar to those given in equation (18) for 
steady flow. They are also derived in Appendix D and are 
expressed as 

4> - K 4> = - K 0” - iv<f (19) 

X X 


where 


y+i 


K* = K + ii^ 





"^Ix "‘^Ix J 


It can be seen that as k or u goes to zero, the oscillatory 
condition reduces to the steady flow condition. This property 
is not present in Landahl's boundary condition (equation 10,11, 
ref, 12); however, for 





^Ix« 


2 

y-1 


and 


k-1,0 

the above reduces to Landahl's equation. These conditions are 
rather restrictive and are not suitable for cases which involve 
flutter or d 5 mamic response analysis , 
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The shock doublet is used in the unsteady case in the same 
manner as is used in the steady case. The same rationale for its 
inclusion is applicable to unsteady flow. 


The form of the shock doublet strength is similar to that of 
a lifting line. In order to be compatible with the pressure 
functions, it is expressed as 


^•-SHq(2) = 


4pU^ 

^AVGq 


SOqi(5) [boqUo(?) 


+ + ...] ( 20 ) 

where / ( 22 ) is the spanwise weighting function given in 
equation (9). The unknown coefficients, bn, are real for 
steady flow and complex for unsteady flow. 


Solution Process with Imbedded Shocks 

For the solution to be rigorous, the upper and lower surface 
solutions for mixed transonic flow should be coupled together 
along the edge boundaries of the wing. Such an undertaking was 
felt to be too ambitious for the current study, hence, the assump- 
tion of decoupled upper and lower surface solutions was retained. 
This assumption was justified through examination of experimental 
data which indicated that the trailing edge pressure coefficient 
does approach zero for small perturbations and it is not too large 
for large perturbations. Thus, the trailing edge and tip pressure 
coefficients were set to zero as is usually done in lifting sur- 
face methods. The leading edge was left alone which resulted in 
some problems in the prediction of large perturbations in steady 
flow. The adequacy of the decoupled assumption will be discussed 
further in the next section on application of the method. 

The potentials, and in equations (18) or (19), are 
calculated as shown in Appendix E for both the lifting surface 
distributions and the shock doublets. Since the change across 
the shock of potential due to interfering surfaces is small, 
its effect is neglected. Hence, <f>'^ and <t>" are calculated only 
from the two surfaces adjacent (upstream and downstream) to the 
shock. The potential derivatives are calculated with the pressure 
coefficient. As shown in Appendix E, the boundary conditions are 
satisfied in steady flow according to the expression 
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Cp(xs,y) - I Gi)(y,y) + //C‘(xs,y) = 0 

and in unsteady flow with 

Cp(xs,y) - Y- GD(y,y) + /^Cp(xs,y) + X G (y,y) = 0 

The Cp (xs,y) and Cp(xs,y) are the pressure coefficient values 
just aft of and forward of the shock. The G"(x,y) and GD(y,y) 
are the chordwise integrals of the supersonic and shock doublet 
products with the kernel function along the downwash chord at 
span station, y. These are the same integrals given in 
Appendix B, 

The downwash induced by the shock doublets is calculated in 
the same manner as is done for surface lift distributions. The 
equations given in Appendices B and C are directly applicable 
where the chordwise integration is simplified to a one-point 
evaluation along the shock doublet. The spanwise integrals are 
identical. 

The physical problem for a wing in mixed transonic flow is 
represented with a "transonic pair" of lifting surfaces as shown 
in figure 8. The upstream surface may have all supersonic flow 
or a subsonic leading edge with accelerating flow to supersonic 
at the trailing edge just forward of the shock. The downstream 
surface must have a leading edge Mach number of less than 1.0 
as well as the remainder of the surface. The shock doublet is 
located between the upstream and downstream surfaces. The 
potential is calculated only from the shock doublet since the 
potential at the leading edge of the subsonic surface is always 
zero for linear theory. 

The aerod 5 mamic matrix construction for a "transonic pair" 
is shown in figure 9. The downwash submatrices An, Ai2, A£i 
and A22 are calculated in the manner discussed in Appendices B 
and C. The matrix Aid is the downwash induced on the supersonic 
surface by the shock doublet and is all zeros as is the matrix 
A12. The A2D matrix is finite and is the downwash induced on the 
subsonic surface by the shock doublet. This matrix is calculated 
as discussed above. The Bsi and Bs2 matrices are the contribu- 
tions due to the supersonic and subsonic surfaces, respectively, 
to the shock boundary conditions. The BsD matrix is the -K<^+ 
term in the boundary conditions due to the shock doublet. The a^ 
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k^‘) 

FOR SUPERSONIC 
SURFACE 


FOR SUBSONIC 
SURFACE 


(-W) 

FOR SHOCK 
DOUBLET 


DOWNWASH ON 
SUBSONIC SURFACE 
DUE TO SHOCK 
DOUBLET 


AT EACH 
DOWNWASH CHORD 


Matrix Construction for a "Transonic Pair" 
of Lifting Surfaces 
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End E2 vectors Ere the usueI pressure series coefficients for 
both surfEces which Ere unknown. The bp vector is the shock 
doublet strength series coefficients. The Wi End W2 vectors 
Ere the known downwEsh boundEry conditions on the two surfEces. 
The iDEtrix equEtion es shown in the figure is solved to obtEin 
the full vector of unknown coefficients, Ei, e 2 End bj). If 
surfEces 1 End 2 Ere not a trEnsonic pEir, then the shock doublet 
End boundEry conditions Ere removed which leEves the usueI inter- 
ference mEtrix form 


\ 

> 

Ai2 


^1 


'^l 

CM 

< 

1 

^22. 


^2. 


W2 


The effect of using the two different solutions on the sEme 
problem will be discussed in the next section on EpplicEtion 
of the method. 


APPLICATION TO TRANSONIC PROBLEMS 


The computer progrEm was Epplied to a vEriety of problems 
involving mixed trEnsonic flow over stEtionEry End oscillEting 
lifting surfEces. The steEdy flow EpplicEtion demons trEted the 
cEpEbility of the method for predicting incrementEl chEnges to 
the SteEdy meEn flow lift distribution. The unsteEdy flow 
EpplicEtions included rectEngulEr End swept trEpezoidel plEnforms 
oscillEting in both elEstic End rigid body motions. The frequency 
vEriEtion rEnged from neEr steEdy flow to moderEtely unsteEdy flow. 
CompErisons were mEde between theory End experiment for uniform 
flow, trEnsonic flow without shocks End trEnsonic flow with shocks. 
In Eddition, a case was considered in which theoreticEl well 
interference effects were Eccounted for which led to En improve- 
ment in Egreement with experiment. Where shocks were used in the 
solution, shock movement wes shown es 2% of the IoceI chord in 
order to illustrEte the resulting lift increment. 
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Rectangular V7ing in Steady Flow 


The transonic method was applied to a case for an aspect 
ratio 3.0 rectangular wing in steady flow at =1.0 (ref. 13). 

The lift change was predicted for an increment in angle of attack, 
Aa=a2 - ai where a2 = 5° and ai = 0°. Predictions were com- 
pared with the sum of the change in lift on the upper and lower 
surfaces obtained from experiment as 


where 


P4a 

L 


II 

n 

- AC 


= 

P«2 

^Pa 


= cL 

- cL 

P4 

Pa 2 

Pa; 


U 

P^a 


da - ct2 “ “i 


Such a definition of incremental lift eliminated error due to non 
zero lift at a=0 and was consistent with a piecewise linear 
representation of transonic flow problems. 


Comparison of prediction and experiment is shown in figure 10. 
Th® fii^st solutxon shown (solid line) is obtained by using the 
average of the predicted pressure changes due to da=5° for the 
local Mach number distributions, Ml, at both a=0° and a=5°. 

For a=0°, the incremental pressures were the same on the upper 
and lower surfaces since Ml was about the same on each (symmetric 
biconvex airfoil). At a=5°, separate solutions were run on the 
upper and lower surfaces. The need for three solutions stemmed 
from the large variation in ML at a=0° and at a=5°. At a=0°. 

Ml varied at midspan from about 0.85 at the leading edge to 1.38 
at the trailing edge. At a=5°, the variation was from about 1.28 
to 1.39 on the upper surface and from about 0.68 to 1.30 on the 
lower surface. In all cases, the control points used were 4 chord- 
wise and 6 spanwise. The problem run time was about 12 seconds 
of CPU time on the CDC 6600 or a total cost of 36 seconds. 
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M«- 1.0 

koo- 0.0 

4a - 5°: a^ - 0°, ai - 5° 

O MEASURED DATA (REF. 13) 

TRANSONIC SOLUTION WITH 

EXPERIMENTAL M). DIST. 
UNIFORM FLOW SOLUTION 


Figure 10, - Results for an AR 3,0 Rectangular 
Wing in Steady Flow at * 1.0 
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A second solution is shown (dashed line) obtained with the 
same computer program for a uniform Ml=0.98 distribution. Com- 
parison of the two solutions with experiment clearly demonstrates 
the superiority of the transonic solution. 


Shown in figure 11 is a comparison of prediction and experi- 
ment for the same configuration as considered in figure 10 but for 
Moo =0.9. In this case a shock existed at about 707o chord on the 
upper and lower surfaces at ct=0O. Again, three solutions were 

for the transonic case in the same manner described previously. 
The control points used were 4 chordwise and 6 spanwise in the 
supersonic regionand 3 chordwise and 6 spanwise in the subsonic 
^®§ion« The total number of unknowns was 42 surface pressure 
function coefficients and 6 shock doublet function coefficients. 

The CPU was 20 seconds on the CDC 6600 or a total of 60 seconds 
for the three solutions. 

Again, the transonic solution is clearly superior to the 
uniform flow solution. Not only that, but the results agree 
remarkably well with experiment. 


Rectangular Wing Oscillating in Bending 

The aspect ratio 3 rectangular wing was treated in the same 
manner as in reference 6 with only a change in the aerodynamic 
method. The Mach number distribution and control point arrays 
are shown in figure 12 along with the results for the (ref. 13) 
case of a=0°, M=0.9, and k^o =0,13 (based on semi-chord). The 
solutions shown are for the transonic method with and without a 
shock and a uniform flow solution. The Mach mamber distribution 
as shown was assumed to be constant in the spanwise direction. 

I^ this case, the presence of the shock oscillation is clearly 
evident in the real part of the solution. Comparison with the 
transonic solution in reference 6 shows a significant improvement 
hhe real part forward of the shock and in the imaginary part 
the shock. The use of shock boundary conditions and down- 
wash point Mach number in the kernel function resulted in correct- 
ing the discrepancies noted in reference 6. 
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Figure 11, - Results for an AR 3.0 Rectangular 
Wing in Steady Flow at Mo© * 0,90 


37 




NC 

NS 

SUPERSONIC 

4 

6 

REGION 



SUBSONIC 

3 

6 

REGION 




CONTROL POINT ARRAYS 




Moo = 0.90 

k 00 = 0.13 (BASED ON SEMI-CHORD) 

OL = 0.0 

© MEASURED DATA 

□ REPEAT 

TRANSONIC SOLUTION 

WITH SHOCK 

TRANSONIC SOLUTION 

WITHOUT SHOCK 
UNIFORM FLOW SOLUTION 


(REF. 13) 


Ml DISTRIBUTION USED 


Figure 12. - Results for an AR 3.0 Rectangular 
Wing Oscillating in M^o “ 0.90 Flow 
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Swept Trapezoidal Wing With Roll Excitation 


This case consisted of a swept trapezoidal wing oscillating 
in a roll mode at 100 Hz in an experimental study conducted by 
Becker (ref, 14), The actual mode included some bending and 
twisting motion which was somewhat difficult to determine from 
the data available in reference 14, This problem was discussed 
in reference 6 due to the fact that the mode shape used would not 
produce a uniform flow solution that would agree with a solution 
given by Becker, As a result, a study was made to determine the 
mode that would produce a solution that would show satisfactory 
agreement with Becker's data at Mo© =0,8, These results are shown 
in figure 13 with a tabulation of the mode used. The slopes were 
all assumed constant in the chordwise direction which probably has 
some influence on the results to be given. 

The first transonic case for =0,937 and k*, =0,61 (based on 
semi-span) was solved with the downwash point array and shock 
geometry shown in figure 14, The results are shown in figure 15 
where the solid line is the transonic solution with a shock and a 
dashed line is the uniform flow solution. The solution given in 
reference 6 agreed better with experiment aft of the shock for 
both real and imaginary parts; however, those results were ob- 
tained with the questionable mode shape as was illustrated by a 
gross overprediction of the imaginary part forward of the shock. 

As was discussed in reference 6, the experimental data were felt 
to be strongly affected by wind tunnel wall interference. These 
effects would be stronger aft of the shock and weaker forward of 
the shock. Thus, the improved agreement forward of the shock 
indicates that the mode shown in figure 13 is more correct. If 
the Mach number increases so as to move the shock further aft and 
increase its strength, the interference effects should decrease. 

The next case for M©© =0,997 and k©© =0,58 was solved with the 
downwash point array and shock geometry shown in figure 14, The 
results in figure 16 show a significant improvement over those of 
reference 6 which supports the suspicion that wall interference 
effects are present in the data. These results also exhibit a 
reasonable progression for the change in Mach number from 0.8 
to 0,937 to 0.997, 
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EMPIRICALLY DETERMINED MODE 



Moo “ 0.80 

koc “ 0.71 (BASED ON SEMI-SPAN) 

a = 0 

O MEASURED, (REF. 14) 

UNIFORM FLOW SOLUTION 

WITH NEW MODE SHAPE 
SOLUTION FROM (REF. 14) 


Figure 13» “ Uniform Flow Solution &nd Mode Shepe for a. 

Trapezoidal Wing Oscillating in “0.8 Flow 
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Figure 14, - Shock and Planform Geometry for the 
Trapezoidal Wing 
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□ 




Ml distribution used 


Figure 15. - Results for the Oscillating 

Trapezoidal Wing at M«) * 0.937 
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Moo “ 0.997 

koo = 0.58 (BASED ON SEMI-SPAN) 

a “ 0.0 

o LOWER SURFACE T MEASURED 
□ UPPER SURFACE J (REFS. 14,15) 

TRANSONIC SOLUTION WITH SHOCK 

UNIFORM FLOW SOLUTION 


0 , 1.0 

C 


Ml DISTRIBUTION USED 


Figure 16. - Results for the Oscillating 

Trapezoidal Wing at M«. “ 0.997 
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Swept Trapezoidal Wing With an Oscillating Aileron 

A low aspect ratio swept trapezoidal wing with an oscillating 
inboard aileron in mixed flow was considered as the next case in 
reference 6. The conditions were Moo *0.942 and koo *0.591 (based 
on semi-span) and the experimental data was obtained by Bergh, 

Tij deman, and Zwaan (refs. 14 and 15). The configuration and 
control point arrays are shown in figure 17. In this case, it 
was felt that wall interference effects would not be so signifi- 
cant. 


The solutions shown in figure 18 include those obtained from 
the transonic method with and without shocks and uniform flow 
theory. The uniform flow solution agrees quite well with experi- 
ment for points on the aileron. Over the remainder of the wing, 
however, both transonic solutions are far superior. The most 
notable point is the large negative value in the imaginary part 
of the measured data at ^*0.55 which is predicted by the theory. 
Comparison with the results of reference 6 shows that the improved 
method is more realistic at all span stations. Inclusion of the 
shock boundary conditions has the greatest effect outboard of the 
aileron. 


Wall Interference Effects in Oscillatory Flow 

The aspect ratio 3 rectangular wing of reference 13 for 
Moo =1-0 ^00 =0.12 was chosen to determine if the transonic 

theory could be used to predict wall interference effects. 

These results are shown in figure 19 at mid-span. The first set 
of solutions shows comparison with the AFFDL sonic box method 
(ref. 1), the uniform flow kernel function solution at Moo *0.98, 
and experiment. With exception of the real part near the leading 
edge, the solutions essentially agree. The next set shows how 
the theory changes when non-uniform flow is accounted for. 

(The Mach number varies from 0.80 at the leading edge to 1.39 
at the trailing edge.) The last set shows the effect of includ- 
ing the first reflection images above and below the wing for non- 
uniform flow (asstiming 1007o reflection). 

The solution is most improved over the aft two- thirds when 
non-uniform flow effects are included and over the forward third 
when wall interef erence effects are added. Since the forward 
third corresponds roughly to the subsonic region, it is expected 
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Figure 17. - Shock and Planform Geometry 

For a Trapezoidal Wing With an Aileron 


45 




Figure 18. - Results for a Trapezoidal Wing With an 
Oscillating Aileron in 0,942 Flow 
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that in cases where shocks are present, the subsonic region 
aft of the shock would be greatly influenced by interference. 
This is felt to be the case for the experimental data given in 
figure 15. It also appears that the transonic results in 
figure 12 would be improved by including wall interference. 


CONCLUSION 


An improved kernel function method has been presented in 
this report that is applicable to multiple lifting surface 
problems in subsonic, transonic and supersonic flow. The method 
yields solutions for steady and oscillatory perturbations about 
linear or non-linear steady mean flows. Results have shown that 
the use of the local Mach number to calculate the kernel function 
leads to significantly improved solutions. The oscillatory normal 
shock boundary conditions as derived have shown the correct limit- 
ing properties as the frequency approaches zero or infinity. For 
cases with shocks imbedded in the flow, inclusion of the shock 
boundary conditions also led to improved solutions. Inclusion of 
an image to approximate wind tunnel wall interference was found 
to improve the agreement with experiment for a case involving a 
rectangular wing oscillating in Moc “1.0 flow. The comparison of 
transonic solutions with linear theory results and experiment 
indicated that the effect of using transonic aerod}mamics in 
flutter and d 5 mamic response analyses could be significant. 


General D 5 mamics Corporation 
P. 0. Box 748 

Fort Worth, Texas 76101, September 30, 1976 
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APPENDIX A 


THE NONPLANAR KERNEL FUNCTION 


The nonplanar kernel function is expressed as 


K(Xo,yo,Zo»^»M) 


g-ikxo r K2T2 

" ’TT' » 


9 


R^^O 
R^< 0 


where (see figure 3 for definitions of and 0q) 

Xo = x-4 

Xo “ y-'n 

2o “ 2 -^ 

“ Vo^ + Zo^ 

R^ = Xq^ + 

= 1 - m2 

T^ = cos(dp-^q) 

^2 “ (Zocos0p-yosin^p)(zocosdq-yosin^q) 

which is valid for either subsonic or supersonic flow. 

The k and M values are those at the downwash point. 

The distinction between subsonic and supersonic flow is 
embodied in the Ki and K 2 terms. 

For supersonic flow (or supersonic downwash points), 
the term is defined as 


(Ala) 

(Alb) 


Ki - Kii + Ki2 


(A2) 
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where 


V /Xo . “ikrui 

Kll = (-^ + 1) e 1 


- 1) . I, 


R - (x<2+^V)^ 


ui - 


Xn - MR 


Xq + MR 
(P-r: 


nc+ikr 


-ikru]^ 


Ef 4 

nc-ikr 


[•- 


“ikrui o ^ ~| 

- 1 + (kr) 2^ („c)2+(kr)2 J 

n=l 


= ikre 


11 . n 
-ikru2 ^^2 

jL^ nc+ikr 
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where 


El = l“l' 


E 2 = e 


-CUr 


c = 0.372 


For the K 2 term, the following expressions are used for 
supersonic flow: 


K 2 = K21 + K22 


where 


K 21 - 


2(il2 + 1) + i£ 

^ ikrM2 


1 R r2 

R (l+u2)^ 



K22 - - 3 Ii 2 - I22 ■ 

For U]^ > 0 , 

u 2n 

-ikrui Pn^l 

I 21 - ikre / ^ 2nc+ikr 

n=l 


2(il£ - 1) + 

)S2xo . 

ikrM^ 1 

R r2 

R 

J 


or for uj^ < 0, 


I 21 = ikre 


11 ir2n 

-ikrui 

2nc-ikr 

. n»l 


+ 2 


- 1 + 


11 




h 

2J./^,v.^2 


(2nc)-^+(kr) 


(A4) 

ikru]^ 

(A4a) 

-ikru2 

(A4b) 


(A5a) 


(A5b) 
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and 


2n 

lo , - lkre '^’^’^2 bn ^2 

^ 2nc+ikr 


n-1 


(A5c) 


The and coefficients in the series summations in 
equations (A3) and (A5) are given in table Al, The a^^ set 
are those given originally by Laschka (ref. 17) for the 
approximation 


u 

(l+u2)% 



-ncu 

an e 


The bn set are those given by Cunningham (ref. 6) for the 
approximation 


(1+u^)^^^ 


1 « 


11 


E 

T-l =1 


bn 


„-2ncu 

e 


The expressions for Ki and K2 in subsonic flow (or at sub- 
sonic downwash points) are similar to but simpler than those for 
supersonic flow. The terms are 


Ki = Kii 

and 

K2 = K21 

where Kn and K 21 are defined for supersonic flow in 
equations (A2a) and (A4a) , respectively. 


(A6) 


(A7) 
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TABLE Al. APPROXIMATION COEFFICIENTS USED IN THE 
KERNEL FUNCTION INTEGRALS Ki AND K2 


n 



1 

0.24186198 

3.509407 

2 

-2.7968027 

-57.17120 

3 

24.991079 

624.7548 

4 

-111.59196 

-3830.151 

5 

271.43549 

14538.51 

6 

-305.75288 

-35718.32 

7 

-41.183630 

57824.14 

8 

545.98537 

-61303.92 

9 

-644.78155 

40969.58 

10 

328.72755 

-15660.04 

11 

-64.279511 

2610.093 
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APPENDIX B 


CALCULATION OF SELF INDUCED DOWNWASH 


Evaluation of the integral given by equation (1) is performed 
in two steps. First the chordwise integration is performed in a 
manner which varies according to whether the Mach n\amber is sub- 
sonic or supersonic. The spanwise integration is then performed 
independently of whether the flow is subsonic or supersonic. 

This appendix summarizes the manner in which these integrations 
are accomplished for equation (1) for q=p, the case of self 
induced downwash, A special treatment is given to chordwise 
integration in the vicinity of the downwash point. 


Chordwise Integration 

For steady supersonic flow, the integral equation may be 
cast in the form 


w(x,y) 


-b Imc 

-2/ / ^C„(|,s) 

5a -1 *’ R So(2-J)2 


where 

w(x,y) 

«2 


Stt w(x,y) 

U 

+ /S^r^ 


(Bl) 


Imc “ Mach cone boundary in the 

transformed coordinate system 

Since the kernel function is singular along the Mach cone 
boundary, £ ^mc» (E )“^ for E-^0, a new chordwise 
variable of integration, 6 , will be defined as 


and 


6 
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1+imc 

Mi 

^■*"imc 
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where 


^mc” Mach cone boundary in the physical coordinate system. 

Substituting 6 for £ and multiplying and dividing the integrand 
by (1- equation (Bl) becomes 


w(x,y) = 




5b 


(^■*"lmc) 


la. 


f 

-1 


2 

( 1 -n 


a-tb 


^p({,5)l 


Xodd 


and 


b(5)d2 

(Y-S)2 


Imc^ 1 


% 


w(x,y) = “ 1“ / / (1-i^) ^ 

° 5a * -1 


(l-i2)^4Cp(|,5) 




b(5)d3 

(y..5)2 


Imc ^ 1 


(B2a) 


(B2b) 


It should be noted that regardless of whether the Mach cone is 
forward ( ^mc<i) (£mc^l) of the trailing edge, the 

chordwise weighting function, (l-‘ 6^)"^ or (1- £2)“%, is the 
same. As a result, the same chordwise integration scheme can 
be applied. 


The Tschebychev-Gaussian quadrature integration formula, 
which is applicable to the chordwise integral in equations (B2) 
is 

1 J 

/ f«)di = 7 f(fj) 

-1 •’ j=i 

where 




and f ( ^ ) is expressible as a polynomial of (2j-l) degree or less, 
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Now, the chordwise integrals in equations (B2) may be evaluated 
in steady or unsteady flow as 

J 


(y-V) K(x-^j,y-77,0,k,M), -1 < 1 


(B3a) 


where 


ij “ - cos(-^^Tr), j=l,2,...j 
ij = % I (l+|mc)^j + (|mc-l) j 


and as 




(y-^)^K(x-^j,y-?7,0,k,M), imc 1 


where 


(B3b) 


ij 




J 


57 



The value for J varies according to the downwash point location 
as follows : 


J = 4+i, i = l,2,..,mq 

where i is the chordwise downwash point number starting with 
i=l at the leading_edge and i^q at the trailing edge. Thus, 
as an example for mq=5, J varies as 

J = 4+1 = 5 (leading edge downwash point) 

J = 4+mq = 9 (trailing edge downwash point) 

The same value of J is used in equations (B3) at all integration 
chords with exception of those near the downwash point for any- 
given downwash point. When the >fech cone falls forward of the 
leading edge, the value of G(^ , v) is set to zero as shown in 
equation (B3a) , ” ^ 


The chordwise integration for subsonic flow is similar to 
that presented above. The integral is of the form 


w(x,y) 


1 


/ 




b(J7)d|d2 

So(r?)2 


(B4) 


The most notable difference between equations (Bl) and (B4) is 
the limits of integration. The chordwise integration formula 
is also changed to 





j=i 


a-ij)f(ij) 


(B5) 


= - 


cos(|^7t), j.=l,2,...J 


The integration is performed by multiplying and dividing by 
(1- f )^/(l+f )^ and using equation (B5). The result for steady 
or unsteady flow is 
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J 

j=l 

(y-»7)^K(x-^j ,y-;?,0,k,M) (B6) 

The value of J is constant at 
J=mq 

for all integration chords except those near the downwash point. 


Spanwise Integration 


In supersonic flow, 


5b 


the spanwise integral is of the form 



(B7) 


where the limits of integration, 1 and 5 ^* the locations 
of the intersection of the Mach fore cone with the wing leading 
edge. Except for a weighting function, (1- 5 ^)%^ and the limits 
2a and 2b> equation (B7) is identical to the equation for 
subsonic flow developed in reference 7. Since equation (B7) is 
more general, it is the form used in this report for subsonic 
as well as supersonic flow. 


Let 


G(y,2) * G(y,y) + (2-y)G'(y,y) + (2-y)^G" (y,y) 

+ . . . 


(B8) 
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Adding and subtracting equation (B8) from (B7) yields 


Ky) 



G(y,?) - (G(y,y) + (^-y)G'(y,y) 


+ (5-y)^G"(y,y) + 

5b 


dg 

(Z-5)2 


+ 


/ 

5a 


G(y,y) + (5-y)G'(y,y) + ( 2 ”y)^G"(y,y) + ... 


dg 

(y-5)^ 


(B9) 


where the first integral may be evaluated with numerical quadra- 
ture integration techniques and the second integral evaluated 
analytically. 


The quadrature integration is performed with a formula 
derived from the equations for function approximation with 
Tschebychev polynomials of the first kind. The approximation 
to a function g( g') over the interval (-1 < g' < 1) is 


S'-l 

g(^') = + X] Ck Tk(g'). 

k=l 


(BIO) 


where 


S' 

Ck = fi f(5;> Tk(g^) 

S=1 


g; = - cos(l|^7r)' , s=l,2,...S 


(Bll) 


S' = total integration chords in the 
range -1 < g < 1 
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If S( 3 ') Is expressible as a pol 5 momial of degree (S'-l) or less 
over tfhe interval (-1 < 2'^ 1), the approximation is exact. The 
quadrature integration formula is derived by integrating equa- 
tion (BIO) from 3 '= -1 to V' = +1 and inserting equation (Bll). 
The result is 


1 S' 

/ g(r’)<i 2 ' • F ^ 8(5^) h (jp 

where 


S'-l 

h(3p =1 S 

® k "0 

and 

= 1.0 


^ = 0 , k*l,3,5, . . . 

1 

“ S Tk^^')^2* “ (k-l)(k+l) > k=2,4,.,. 
-1 


(B12) 


The only difference in the quadrature integration formula above 
and that used in the subsonic spanwise integral of reference 6 
is the definition of h(]?g). 


In order to achieve a more optimiam distribution of spanwise 
integration points, a coordinate transformation is made such that 
the greatest density of integration points is next to the down- 
wash chord. This step is particularly necessary for downwash 
points near the leading edge in supersonic flow. The transforma- 
tion is 

!?s “ Yr + (5s+l) » ?s < “ yr 
2?s = Ir + (!?s-l)> yr 

where 

yr = cos(ll|?-) , r=l,2,...n , n=l,3,5,... (B13) 


61 



and 


S = S* = n(2n+l) type 1 wing in figure 2 (Bl4a) 

S = 2S' = 2n(TT+l) type 2 wing in figure 2 (Bl4b) 


S - S' = n(n+l) type 3 wing in figure 2 (Bl4c) 

The yr defined in equation (B13) are the same as used in subsonic 
flow~and__are interdigitated with the Vs points by the relationship 
between n and S in equations (Bl4). 


Returning to equation (B9) , the integral is evaluated with 
equation (B12) as 


S 


I(y,) - f L HVi) + Q(^r) 

® s-1 (yr-3s)^ 

where 


(B15) 


Q(Yr) = G(yr,yr) Eofe) + G'(Yr>Ir) 


+ G'*(yr>yr) E 2 (yr) 


(B16) 


The En(yr) terms are the differences between the quadrature and 
analytic integral evaluations 


Eo(lr) 


Ei(yr) 


1 

. (2?a-Yr) 


1 

(^b-yr) . 


S 

S (?s-Yr)^ 




2?b-yr 

!?a-yr 


7T 

s 


s 

V h(!?;)d(5g) 

(Vs-yr) 


s 

E2(Zr) = (!?b-2a) “ 4* h(v^)d(Vs) 

^ s=l 


(Bl7a) 


(Bl7b) 


(Bl7c) 
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for 


^(2s) “ 1 » 5a^2?s^5b 

* 0 , Vs<V& or 5s>2?b (Bl7d) 

where in the E© term the Mangier formula was used (ref. 18). 

No terms higher than E2 (^r) have been found necessary, hence, 
the series was truncated at that point. 

The same equations are used in subsonic flow by simply 
setting the limits at 

3a = -1 

3b ■ +1 

which maintains continuity from subsonic flow to supersonic flow. 

Shown in figure Bl are examples of the relationships between 
7r» ’?s > !?a !?b* Also shown are hypothetical distributions 

of G(yr> ^ )/(yr"t?s)^ which are being integrated. The two 
examples illustrate the differences between the treatment of 
type 1, type 2 and type 3 planforms shown in figure B2. The 
differences are embodied in the definition of the spanwise 
functions and v : 


V = 2So 


type 1 surface 


(Bl8a) 




f£ 

2 


(1+5) 


type 2 and 3 surfaces 


(Bl8b) 


^ (1+2) b 3 rpe 2 surface (B18c) 

For equation (Bl8a) , there are S* total integration chords on the 
wing and n ■ (s'/n-l)/2 downwash chords on the right hand wing. 
For equations (Bl8b) and (Bl8c), there are 2S' and S' total inte- 
gration chords for type 2 and 3 surfaces, respectively, and 
n “ (S* /n-1) downwash chords on the right hand wing. 


The G, G' and G" terms are simply calculated as the 
coefficients of a quadratic curve fit of the chordwise integrals 
at G(yr»3?S')» G(vr»^r)» G(vr, ’?!)• The relationship between 
these integrals is also shown in figure Bl. In unsteady flow, 
these terms will have real and imaginary parts. 
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Figure Bl. - Geometric Relationships in the 
Spanwise Integration Scheme 
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Chordwise Integration Near and 
Including the Downwash Chord 

Because of the close spacing of the adjacent integration 
chords with the downwash chord in the spanwise integration scheme, 
the chordwise integrations must be extremely accurate in this 
region. Thus, all chordwise integrations within a small spanwise 
distance of the downwash chord, yr, are evaluated more precisely 
than for those chords outside of'this distance. The distance is 
taken as +(l-yi) where yi is defined by equation (Bl3). Thus, 
for all Tjg defined as 

(jr“4) <Us< (yr+4) 


where 


A = 1-yi 

the nximber of chordwise integration points, J, are defined as 
j' = 3J+1 subsonic integration 

T I J- O T 

supersonic integration 

The terms G(yr>I?s)» G(yr>yr) G(yr» ’J l) always included 
within this region which results in very accurate estimates of 
the G, G* and G” terms required for equation (B16) for subsonic 
flow. For supersonic flow some additional manipulation is needed. 

Because the chordwise integrand is singular at ^ = ^tnc 
for ^ 9 ^ yr, and is non-singular f or 5 = yj., a difficulty arises 
in the calculation of a value for G(yr,yr) that is continuous 
with G(yr» 2 J ) as As a result, a simple but devious means 

for calculating G(yr»yr) as a function of G(yr» V ) '^as developed 
as follows. For simplicity, let 

Gg “ G(yr»!?s) ^1 “ 

Gr ” G(^r,yr) 

As a first approximation, Gro» to the value of Gr without actually 
calculating the integral, an average value is used 



(B19) 


65 



Since the streainwise variation of the kernel function becomes 
rapidly independent of ^ for 5 w yr, the Gr© approximation 
actually represents a chordwise integral of the form 

/ ^Cp(|.yr)(K-f) 

^LEr 

where 

= x-^d' 

= yr“"?S = 

and f T,F-r is the leading edge location at 


+ (1-m2)(4')^ 




d^ 


(B20) 


The exact value required for Gr is 

X 

/ ACp(f,7r)df 
LE 

Thus, a correction must be made to Gro* The first correction is 
made by removing the singularity portion at ^mc without disturb- 
ing the leading edge sweep effect as given by ^lE. Since the 
error in equation (B20) is near the aft limit, , the next 

level of approximation, Gri, is ^ 1 


^^1 + ^^ACp(|4',yr) 

- •[ (x-0 
^LEs 


/ 

^LEr 


d^ 


(x-0^ + (1 -m2)(^')^ 




11 

2 


/ 


(X-^) 


’LEi 


(x-0^ + (1-M^)(A’)^ 


1 “% 


if 

2 


which becomes 

Gri “ Gr^ + ,yr) 


(^4-^LEr) “ 7 




+ (1-M^)(d')' 


% 


1 , 

2 




k 
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Now, Gri is approximately 


Gri - f ACp(f,yr)d^ 

^LE 

Completing the integral from to x yields the final approxi- 

mation to Gr 


Gs-H3i 2b (2) _ 

Gr -^5— + -^-" 


2 Sr, ACpCf^'jJr) 




1 

2 


.2 , /T t^^2\ / a! \ 2 


(’=-^LEs^ + (1-M^) (/!')■ 


1 

2 


(==-^LEi)^ + (1-m2)(4')2 


+ 


b(:?)y3i* 


ACp( jXr) "^^pC^jyr) 


(B21) 


The use of equation (B21) in place of (B3a) to obtain Gr in 
supersonic flow results in the desired continuity in the chord- 
wise integrals used in the spanwise correction terms. 
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APPENDIX C 


CALCULATION OF INTERFERENCE EFFECTS 


The equations presented in Appendix B are used to evaluate 
the self-induced downwash in equation (1) when q=p. This Appendix 
summarizes the equations for calculating interference downwash in 
equation (1) when q?^p for which the procedure is essentially the 
same as given in Appendix B. The chordwise integration is per- 
formed first which for supersonic flow requires a special treat- 
ment for a 3/2 power singularity in the integrand along the Mach 
hyperbola. The spanwise integral follows which requires treament 
of the non-planar kernel function such that the integral converges 
to the coplanar case as (z 0. 


Chordwise Integration 


The form of the integral equation for parallel non-planar 
surfaces in steady supersonic flow is 

^mc 


Wp(x,)0 

u 


47TPU^ 


/ 


/ 

Ile 


2x, 


ACp(|,2) 


r^R 


1 - 


2z 


2 „2i 


zTjS 


_ o 


d^drf 


(Cl) 


where 77 a and are the left and right hand limits of integra- 
tion defined by wing geometry and the Mach hyperbola. If the 
chordwise integration is carried out to the Mach hyperbola, 
a 3/2 singularity is encountered of the form 


LIM 

^->^mc 


(x-^)2-ys2r2 


1 -3/2 


Hence, the finite value of the improper integral must be taken. 

In order to determine the finite value of the integral, 
the differentiation of the potential equation for the downwash 
must be considered. Let the downwash be expressed as the 
following simplified form: 

’^b ^mc 

^ f f ^Cp(^,T?)<^(x,^,y,T7,2,o)dfdi? 

^a ^LE 
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where 0( ) is the potential at point (x,y,z) due to a point 

load at ( ^ t rf y 0). Since the differential operator is outside 
the integral sign, Leibnitz's rule must be used to perform the 
operation 

b(t) 

Ie / g(x,t)dx = g(b(t),t| 

a(t) 


b(t) 


-g|a(t).tl / i|g(x,t)|dx 

a(t) 

The limits a(t) and b(t) are constant in subsonic and coplanar 
supersonic flow; thus, their derivatives are zero. For the 
supersonic non-coplanar case, however, the derivatives are not 
zero. Since the chordwise term is of much greater importance 
for interference effects, only the derivatives resulting from 
the variation of ^mc with z will be accounted for. 


For the cherdwise integral, H(y,i 7 ), a constant pressure 
distribution, ACp, will be assumed for simplicity. The integral 
may then be expressed as 

^mc 

H(y,T7> = dCp ^ J* ‘^(x,|,y,i?,z,o)d^ 

^LE 


which, with the application of equation (C2) , becomes 


4Cp 




^^mc 

dz 




For 




2U ^0^0 

■ 4:7/) r2R - *oS/3r 


it can be shown (after some algebraic manipulations given in 
reference 6) that the following is true 


^^mc 


2 2 

2u r 


‘^(x,4^mc»y,^>z,o) = - ^ “72 


LIM 


(C4) 


The differentiation inside the integral in equation (C3) is the 
kernel function 


^ t s 2U 5 Cq 


1 - 


2 2 ? 
2zo Zo r 


■ 4^ K(x.«,y,i),2,o) 


(C5) 


singularity yields 


H(y,r,) = 


dCpU 

^7rp 


T' 


r" 
^mc 


+ / 

^LE 


(C4) and (C5) and 

removing 



^mc 

LIM 1 


f ^ 

+ 

J 2 

^LE ^ 

K(x,^,y,’?,2,o) - 

2^0^ oV ' 

o o 


2,, 2 


d^ 


r^R- 




(C6) 


But 


^mc 


d^ 


2zo18 


2«2 


LIM A. 


^LE 
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Thus, equation (C6) becomes 




4CgU 

^TTp 


2zq2^2 

r^RLE 




me 



K(x,|,y,i7,z,o) 


2xoZo2ys2 

t2r3 




which is the final form for the chordwise integral with a 
constant pressure distribution. For a variable ACp( i » -q ) ^ 
the singularity is treated at the Mach hyperbola, hence. 


H(y,^) = 


U 

47rp 




^Cp (|mc j^) 


^mc 


+ f ACp(f,5)K(x,^,y, 17,2,0) - ^XqZo /3 ACp(|mc,!?) 

^LE L ^ 


d^ 


Thus, for steady or unsteady interference of non-planar surfaces 
the quadrature integration of the chordwise integral becomes 

J 


b(22) (l+|mc) 7 t 9 ^ _ 

Hq(Z»2) = So 2 ACpq(|j,‘5)K(x-^j,y-i7,z-^,k,M) 

j=l 


T2 ACpq(fmc,5) 


i , b(r,) i 

RlE ^ “ r.3 


j=l "J 


, ^c ^ 1 


where T2 is defined for equation (Al) and 

2 2 2 2 

Rj = - pV 

= - cos (^^^) , j=l,2,...J 
Ij “ % I (1+lmc) + (fmc"l) | 

(!2) I 
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For the case of the Mach hyperbola falling aft of the trailing 
edge, the 3/2 singularity need not be considered. Hence, the 
chordwise integral takes on the form given in Appendix B 

J 

h 

Hq(y.5) = ^ T ^Cpq(£j,?) . 

j=l 

K(x-^j ,y-Tj,z-C,k,M) , fmc^l (C7) 


where 


= - cos(^|~7t) , j=l,2,...J 


For subsonic flow, the chordwise integral is similar to 
equations (B6) and (C7) 

J 

Hq(Y,!2) ^Cpqdj,^) . 

j=l 


K(x-^j,y-T?,z-4,k,M) 


(C8) 


where J = mq and 


For »? = y, 

J " 3mq+l 

as discussed in Appendix B. 
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Spanwise Integration 


The spanwise integral for non-planar interfering surfaces 
may be written as 

Vh 

= / . [ Tl 2T2 




!?a 




drj 


where the sign / denotes a "Pseudo Mangier" evaluation of the 
integral as discussed in reference 6, The Gq(y, r) ) function is 
the modified chordwise integral 


Gq(y,2?) 


Hq(y,^) 

Ti - 2T2' 

r2 


which will be expanded as 


Gq(y»!?) = Gq(y,y) + (^-y)G^(y,y) + ... 


Following the developments of Appendix B, the spanwise integral 
is evaluated as 

S 

Iq(Zr) = h(5s)Gq(yr,5s) + Qq(Zr) 

s=l 

where h( 17 s) and 77 g defined in Appendix B 

Qq(Zr) = Gq(yr,yr)|Fi(yj.) cos(0p-0q) 

+ F3(^r) &in(0p-9q)| 

+ G^(yr»rr)|F2(Zr) cos(0p-0q) | 
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and 


Fi(^r) 


:^r-Vb _ Yr-la 

r 2 r 2 

rb ra 


S 




d(3s) 


f2(iTr) = •Ia 



+ 2 ' 


r 2 
*■ rb 



S 

“s" ^(!?s) (5s"Zr) 

s=l 



^ 1 
r 4 
rs-* 


d( 2 s) 





TT 

s 


E h(2s) d( 3 ,) 


S =1 


r 4 

rs 


where 


<l(2s) “ 1 » !?a^ 5s^2?b 
<3(!?s) = 0 . 5s'< I/a , !?s< Hb 
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The variables (^, ^ ,z) are all measured relative to the surface 
over which the integration is being performed as shown in 
figure 3. The perpendicular distance from the surface to the 
downwash point is ^ and it is measured at span station y. The 
terms ^rb> ^rs 

2 2 2 
^ra “ (Yr“!?a) + S 

^ 

^rs ~ 0fr"!2s)^ 5^ 

For subsonic flow, 

2a = -1 

2b = +1 

This completes the equations necessary to perform the spanwise 
integration for interference effects due to a strearawise planar 
surface of arbitrary orientation in steady flow. 

For unsteady flow, it is necessary to include additional 
terms which account for the logarithmic singularity in the 
spanwise integrand, Vi^JU\x\ . This singularity is relatively 
weak for most problems, i.e,, k<l. The correction terms 
developed in subsonic flow (ref,. 7) are applicable in this case. 
The spanwise integration function, h( ^g), is given in Appendix B, 
equation (Bl2) . The form for the unsteady correction terms 
becomes 




k2 

+ cos(0p-0q) 


Gq(^r»yr) 


Ji(yr) 


S 


- f- X d(2?s)h(2s)^rjg 

S=1 


+ G^(Zr.yr) 


J2(yr) 


s 

S=1 
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where 


Jl(^r) * (yr-2a) ^ ^ra " (^r'^b) 


+ 2|so| 


J2(yr) - - j r|a^ rjg - r|b^ + ?a+3b-22r 


which completes the equations for unsteady flow* 
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APPENDIX D 


THE NORMAL SHOCK BOUNDARY CONDITIONS 
ON SMALL FLOW PERTURBATIONS 


For a mean flow in which strong shocks are imbedded, special 
boundary conditions must be satisfied across the shocks on small 
perturbations to the mean flow. This Appendix presents the deri- 
vation of the shock boundary conditions first for steady and then 
for unsteady flow perturbations. These boundary conditions 
satisfy the Rankine-Hugoniot condition and continuity of potential 
across the shock. 


Normal Shock Boundary Conditions for Steady Flow 

Referring to figure Dl, the upstream conditions are M", 
and U" (referring to Mach, total potential, and velocity normal 
to the shock) and the downstream conditions are M+, and U+. 

Beginning with the Rankine-Hugoniot relation as given 
by Pai (ref. 19): 


or 


= 2:1 + 2 

U" y+l (M")2(y+1) 


u^ = u' 


-j_4. 2 

>'+1 (M“)2(y+i) 


The velocity change is obtained 


(Dl) 


AU = 


U“-U"'‘ 


2U- 

y+1 


1 - 


(M-)2 


2U~ 

y+1 


c"\2 


-(f) 


(D2) 


Now, let 


U“ = 1 + (/)■ 


From which 


= 1 - 2 


(D3) 
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Dl,- Flow Variables in the Vicinity of an 
Idealized Normal Shock 
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From Landahl (ref. 12, equation (1.4)), can be expressed as 

.2 _ 1 


^ - (y-l)c^^ 


X 


for first order effects. Thus, 

- (y-m 


m2 


X 


y: 


1 ” 2^jj + ... 


from which the first order terms are 



1-|2 + hI 


For 

u+ - 1+4 

the velocity change is 


(D4) 


4u - (1+4) - (1+4) - 4 -4 

Substituting equations (D3), (D4) and (D5) into (D2) yields 


(D5) 


(/>; - c|,+ = -1- 

^ ^ y+i 




1 -i 


m2 


- ( 2+M2(y-l))4 


from which the first order terms are 


- 4 


2 

7+1 


(D6) 


Letting 


4 = 


7-1 

7+1 


(7+1)m| 


equation (D6) becomes 


(D7) 


<ii + Mc/)“ = 


X 


2 

7+1 


(i - 1) 

U2 ' 


(D8) 


which differs from the derivation of Landahl by the presence 
of iJL , For 1^1, M=l> then equation (D8) would correspond to 
Landahl's equation (equation (10.9), ref. 12). 
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The total potentials, 4>~ and 0"**, are the stim of the mean 
potential, plus the perturbed potentials, </>, 

4>~ = (f>l + + <#>’ (D9a) 

= *+ + + ,),+ (D9b) 

where X is the shock displacement as shown in figure Dl. 

From the condition of equality of potential across the shock, 

0" = 0+ (DIO) 

The following relation between X and <f> 


X 


* - ~ 


(Dll) 


since 

<^1 = 

The fundamental assumption which underlies equation (D9) is that 
the potential 0^ is analytically continuous on both sides of the 
shock, and that it can be expanded in a Taylor series forward 
and aft of the shock » Equations (D9) contain only first order 
variation of 0 due to a perturbation 0, therefore, equation (Dll) 
is the first order movement of the shock due to 0 that is 
necessary to maintain equality of potential. 


Another relationship can be obtained between X and <f> by 
differentiating equations (D9) 


^Ix ^*^lxx ^x 


0"*" = + X0^ + ^ 

X lx Ixx X 


and substituting into equation (D8) . Thus, 
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which separates into the following for the zero order terms 


^Ix + 


= —2— / JL. - 1 ) 

lx y+1 ' ]yj2 / 


oc 


and the perturbation terms 


X = - 




Equation (Dl3) is the first order shock movement due to (f> 
that is necessary to satisfy the Rankine-Hugoniot relation. 


(D12) 


(D13) 


With equations (Dll) and (D13) , 
to obtain 


can be eliminated 


lx lx 




which may be simplified to 


where 


- K<^ 

X X 


K 



+ 





(D14) 


(D15) 


Since the constant K is a function of the mean flow potential 
and <f>“ is unaffected by conditions downstream of the shock, 
equation (D14) represents a constraint on the downstream 
perturbation 
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Normal Shock Boundary Conditions for Unsteady Flow 

The Rankine-Hugoniot relation as given in equation (D2) is 
also used for unsteady flow. The velocities, U" and U^, are 
expressed as 



ikt 



ikXe 

(D16a) 

= 1 + - 

, ikt 
ikXe 

(Dl6b) 


which includes the shock velocity in the exponential terms, 
Landahl has given an expression for the speed of sound in 
unsteady flow (equation 1.4, ref, 12) as 


2 

c 



(r-l)Wx + i(>t) 


with equations (D16) and (Dl7) , the expression for (1/M")* 
becomes 


(D17) 


( 




- ikt. 

. ^oc 

1 - - ik\e ) + ... 


from which the first order terms are 
-\2 , 


c^ 

u- 


m2 

lioo 


2/^_i\\ j. ~ t vr2/-- t X i“ .I . ikt 


1 - (2 + M^(y-l)) + 2ik\e 

Substituting equation (Dl8) into (D2) yields 


(D18) 


(C - ^ 


y+1 


- ikt 1 

1 + <t>^ - 


1 - i 



^OC L 


1 - (2 + Mi(y-i)) 0 


X 


- M^(y-l) + 2ik\e 


ikt 
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from which the first order terms are 




^ / 1 . 1 

y+1 y+1 ^m2 " 


where 


r= 


2 (y- 1) 

y+1 


(D19) 


Except for the exponential and terms, equation (Dl9) is 
identical to (D8) . 

Following the development for steady flow, the total 
potentials due to an unsteady perturbation are 

+ (fi (D20a) 

<j>^ ^ <i>^ + + <t> (D20b) 

The equality of potential condition for unsteady flow yields 
equation (Dll) for the shock movement necessary to maintain 
this condition. 

Next, the derivative terms for equation (Dl9) are obtained 


from equation (D20) 





''’ixx + 

ikt 

d> e 

(D21a) 

^ + Xe^^^ 

X '^Ix 

<xx- 

/ eikt 

X 

(D21b) 

</>t = ike^^^ ( X <t> 

lx + 


(D21c) 
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Substitution of equations (D21) into (D19) yields 


/ ,+ , ^^ikt. + . * + ikt. 

"^Ixx + ® ) 


. /j.’ . V ikt,- 

+ /i((|)ix + <;>ixx + 


iktv 

</>x e ) 




lx 


+ <f) - 


4ik 

y+1 


,ikt 


■ m (i? ■ M 


from which the first order terms obtained are 
^Ix ^^Ix “ y+l (m2 ' ^ ) 

and 


X = 


^x ^*^x 




+ H-<i> 


Ixx 


1 V (f>~ _ 4ik 
lx y+]_ 


where 


2(y-i) 

V = Fk = y_i k 

Equation (D22) differs from its counterpart in steady flow, 
equation (Dl3) , by the terms containing xv , 

With equations (Dll) and (D22) , X can be eliminated 
to obtain 

- <f>~ _ ^x + A^^x + 

Vx - ‘'■ix " Mk ''*Ixx ^ ‘"'’lx ■ ^ 

which may be simplified to 

(f>t - - iy</)‘ 


where 


+ 


K* = 


hxx + + ^^^Ix ~ 

4k - 


4ik 

y+1 


(D22) 


(D23) 


(D24) 
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or 


K* = K + iv 


~ y-1 

- ■#’ix J 


(D25) 


which completes the normal shock boundary conditions for 
unsteady flow. 

Comparison of equation (D24) with the boundary conditions 
derived by Landahl (equation 10.11, ref. 12) shows a significant 
discrepancy. However, with the restrictions 

k ~ 1 


and 

~ (‘^Ixx + m4>1xx) 
equation (D24) reduces to 


<h + 


- - - + 


4ik 


+ 2ik(y-l) 

(y+ r+ 1 




X 


or 

. + 4ik</>~*' _ 4ik</> i 

^ ■ (r+ (T+ 1)«-^ - ^ 

- ^ - ^tx)| 

which is identical to Landahl 's equation. Thus, Landahl 's 
boundary conditions are valid only for high frequency 
oscillations and weak shocks as is implied by the restrictions. 
Equation (D24) is therefore more general and is necessary for 
use with a method which will be used in flutter and d}mamic 
response analysis. 
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APPENDIX E 


CALCULATION OF POTENTIALS FOR THE 
SHOCK BOUNDARY CONDITIONS 


The potentials needed for the shock boundary conditions are 
the supersonic potential just forward of the shock, (f>~ , due to 
the upstream surface and the subsonic potential aft of the shock, 

, due to the doublet. Also needed are the potential derivatives, 
<<>x and <l>^. 


The supersonic potential at point x,y,z due to a doublet 
sheet lying in the ^-r) plane is the integral over the sheet: 

^ ’^a?LE 

The nature of the term 

in the integrand is such that 




\j J.2 
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Thus , 



a 


(y, V )dr} 


(E3) 


where G(y, 17 ) has the same definition as given in Appendix B. 

In order to treat the singularity in equation (E3) , the integrand 
is modified as follows : 


_±_ 

Ubr 


hf r2 I (y-’l) 

1 /» t) 

- i; ° 


- G 



( y. 



The first integral can be evaluated numerically and the second 
analytically as 


/■ 


zdr] 


'n„ 


(y 


. „n-i . e.n-1 


y - 


Taking the limit z-*-0 for both integrals leaves 

_ 1 


Ub, 


8 


G (y, y), Va< y < Vh 


(E4) 


The same equation is applicable for the doublet line. 


The potential derivatives are directly related to the 
pressure coefficient forward and aft of the shock. For <f>x, 
the derivative is taken with respect to the aft limit of 
integration, hence 


U 


I Cp (Xj, y) 


(E5) 


where Xs is the shock location. Equation (E5) is applicable to 
either the supersonic or subsonic regions. Since the shock 
boundary conditions are satisfied aft of the shock, the value 
of Cp for the shock doublet is zero. 
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The normal shock boundary conditions are satisfied with the 
following terms in steady flow: 


^ Cp (Xg, y) 

= ~ I G" (y, y) 

- g I G"(y, y) + Gi)(y, y) | 


pressure coefficient at the 
trailing edge of the super- 
sonic region at y, 

chordwise integral directly 
upstream of the shock boundary 
condition point, (xg,y). 

pressure coefficient at the 
leading edge of the subsonic 
region at y. 

total potential just aft of 
the shock due to the supersonic 
region and the shock doublet. 


Substituting these quantities into equation (18) or (Dl4) yields 
“ 4 Cp (Xg, y) + f [ G-(y, y) + Gj, (y, y) | 

= 4 S y) I y) 

or 


cj (xg, y) - j Gd (y, y) + /xC" (xg, y) = 0 

It should be noted that in steady flow, the G’(y,y) term is 
cancelled. 

For unsteady flow, all of the above terms are the same 
except they are complex. Substitution into equation (19) 
or (D24) yields 

1 + K* t 1 

■ 4 ^^s^ 8 I I 

“ 4 y^ ■’'1' (y^ y) 

or 

Gp (Xs> y) - I Gq (y, y) + ^C’ (Xg, y) + ^ G" (y, y) = 0 
For unsteady flow, the G“(y,y) term is retained. 


(E6) 


(E7) 
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APPENDIX F 


CALCULATION OF LOCAL MACH NUMBER 


The local Mach number, Ml, is used extensively in the tran- 
sonic algorithm. As a result, it must be inexpensive and simple 
to calculate. The input data from which the local Mach nimber 
is obtained must be simple and flexible in format. These require- 
ments are satisfied through a two-step procedure which starts with 
input values of Ml given at even spaced points over the planform 
and ends up with a function for Ml in the form of Tschebychev 
polynomials , 


The input data is given at even spaced chordwise intervals, 
A , starting at the leading edge and ending at the trailing 
edge. The chordwise distributions are given also at even spaced 
spanwise intervals, A^ , starting at the inboard tip and ending 
at the outboar^ tip. The intervals are tied to the total number 
of chordwise, m, and spanwise, n, downwash points as follows: 

A^= ^ (chord fraction) 



There are 


(semi -span fraction) 


m + 1 values of Ml input chordwise at 
n + 1 spanwise stations. 


The locations are the sets 


^i 


= ( 0 , 


1 2 

m’ fn" — 


1 ) 


(Fla) 


’’j ■ n^ h’ 


(Fib) 


It is possible to input a single chordwise distribution into the 
program which is used at all span stations if desired. 
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Next, the input Ml are fit with the surface spline procedure 
described in Appendix G that is used for structural mode inter- 
polationo The surface spline is used to calculate Ml 2 values at 
the optimum fitting points for Tschebychev pol)momials of the 
first kind, Tn. 

The Tn pol)niomials are defined as 

To(x) = 1 
Tl(x) = X 
T 2 (x) = 2x^-1 


Tn(x) = 2x Tn-i(x) - Tn_2(^) 
over the interval 

-1 ^ X ^ 1 

From reference 20, an arbitrary function of x can be expressed 
as an N term expansion of Tn(x) as follows: 

N-1 

f(K) - ^ +^C„ T„ (X) 

s (*■2) 

The coefficients, Cn, can be obtained by a simple vector 
multiplication 

N 

' I S ^ '^n <=^1) 

(F3) 

where the xi are roots of the Tn+1(x) pol 3 momial. 


x^ = cos 


i - 1, 2, 


N 


(F4) 
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Thus, the new values of Ml are calculated at 


xi = (xi, K 2 > ... x^+i) (F5a) 

yj = (yi^ vi’ ••• yfi+i) (F5b> 

where the and yj are calculated with equation (F4) for 
N = (ni+1) and (n+1) , respectively. The relationships between 
the ( ^ , V ) in equation (Fl) and (x,y) in equation (F5) are 



Since 0 < f < 1 and 0 < r? < 1. 

Once the Cn are calculated from equation (F4) with the new 
values of Ml 2» the value of Ml at any point on the planform can 
be obtained with equation (F2) , 
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APPENDIX G 


STRUCTURAL MODE INTERPOLATION 


The structural modes are interpolated in the program with a 
surface spline fit (ref. 21) which simulates the deflected shape 
of an infinite plate pinned at the points which are being inter- 
polated from. The scheme is one of the more dependable methods 
that is currently available. It has some disadvantages, however, 
in that it does not extrapolate well under certain conditions and 
it tends to "sag" in cases where there are large spaces between 
fitting point groups. 


The spline equation is 


N 

n*l 

For n=l,2,«.,N where wj is the function value at point (xi,y^) 
and 

^ni “ + (Xn - 


(Gl) 


The points (xn,yn) are the fitting points. The coefficients 
[ai,a 2 ,a 3 ,bi, . , .bj^] are determined by equating the wj to the 
known deflections at points (xj,yj) = (xn,yn) for n=l , 2 , . , .N. 
Three other equations are also satisfied, 

N N N 

n*l n=l n=l 

to give the necessary N+3 equations to solve for the N+3 unknown 
coefficients. With the coefficients known, deflections can then 
be calculated at any point desired within the planform. 

A scaling and coordinate transformation is performed such 
that the (x,y) values never exceed 1.0. A structural surface 
is defined in different ways as shown in figure Dl. For the 
comer points, xi^, xi 2 , yi and x©]^, Xo 2 » Yq* ^ point (x,y) 
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within the structural surface is transformed to ( ^ , t? ) as 




b(y) 


V = 


y-ym 

S 


where 


Xm = 



y - yi 

2S 


^ V + ^02 - 


X 


b(y) = 




- X. 


ii y - y^ / ^02 - ^oi - xj2 + 


2S 


y ,^1+^0 

■^m 2 

s = yp -H yi 
2 


As can be seen in the figure, use of a square structural surface 
simply scales the (x,y) to ( , i? ) which are geometrically similar 

but less than unity in magnitude. The trapezoidal structural 
surface skews the structural points into a square plane which is 
sometimes desired. Skewing is not recommended for delta or other 
pointed tip planforms. 
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APPENDIX H 


PROGRAM UTILIZATION AND 
INPUT DATA ORGANIZATION 


The program is applicable to steady and unsteady flow in all 
Mach number ranges, subsonic, transonic and supersonic. It can 
be used to obtain uniform flow solutions in both subsonic and 
supersonic flow. Non-uniform flow solutions can be obtained for 
all Mach ranges including mixed transonic flow with imbedded 
shocks. In order to use the non-uniform flow capability, the 
Mach number distributions and shock locations are input and must 
be obtained from an outside source such as experiment or a finite 
difference solution. 

Interference effects can be calculated for coplanar, non- 
coplanar, and intersecting planar surfaces. The aerodynamic 
geometry data format permits the input of arbitrarily arrayed 
surfaces which are constrained only to being streamwise. Each 
surface may or may not have an image surface, hence, asjmunetric 
configurations may be constructed. 

The program can be used to calculate generalized forces in 
unsteady flow and pressure distributions in either steady or 
unsteady flow. The solutions obtained are essentially independent 
in computer cost of the number of modes or downwash vectors input. 
Once computed, inverted, and saved on a magnetic tape, the aero- 
d3mamic matrices can be used on subsequent problems for very 
little cost as long as Mach number, reduced frequencies, and 
aerod5mamic geometry remain unchanged. Thus, the method is 
tailored for design applications where the structural mode shapes 
change continually for structure changes and payload variations 
while aerodynamic parameters remain constant. 

Presented in the following sections are the detailed instruc- 
tions for providing input data to the program. The first section 
describes the general arrangement of the data and the second 
describes the specific items in the data. 
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Deck Arrangement 


The major input data classifications are: (1) the library- 

data which contains structural geometry data and the mode shapes 
and (2) the problem data which includes all aerodynamic geometry 
and option data. The program can be used to calculate generalized 
forces or pressure distributions. The generalized forces can be 
obtained only for unsteady flow whereas the pressure distributions 
can be obtained for steady or unsteady flow. 

Library data organization .- The library data organization 
presented below is in the form of card sets, L-N, where L refers 
to library data and N is a sequence number corresponding to a 
specific function as described. The cards sets may be composed 
of one or more cards and their contents are described in the 
following section. All card sets are required if library data 
is used with exception of L-6 and L-7 which are optional as noted. 


Card Set 

Func tion 

L-1 

"LIBR” card 

L-2 

option data 

L-3 

cons tants 

L-4 

structural point locations 

L-5(l) 


• 

. structural surface data for "NSUR" 

• 

total surfaces 

L-5 (NSUR) 


L-6 

zero deflextion points - supplied if 
option L3?^0 on card L-2 

L-7 

mass data - supplied if L2?^0 on card L-2 

L-8 

mode deflection data 
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Problem data organization , - The problem data organization is 
presented below in a format similar to that of the library data. 
The sets denoted as P-N are mandatory for all problem decks 
regardless of whether or not library decks are used. Those 
denoted as P(L)-N are used only when the library is used and 
P(P)-N are used only when no library is used. Specific contents 
of each card set are described in the section following that for 
the library data. 

Card Set Function 

P-1 title data 

P-2 major options 


P-3 general aerodynamic data, integer 

P-4 general aerod 5 mamic data, real 


P-5(l) 

P-6(l) 

P-7(l) 

P(L)-8(1) 

P-9(l) 


integer control data 

► geometric data 

structural surface used 
if "LIBR" data is input 

Mach distribution data 
if ITRANS?^0 on card P-5(l) 


1 


Data for 
aerodynamic 
surface no. 1 


J 


P-5(NSURF) 

P-6(NSURF) 

P-7(NSURF) 

P(L)-8(NSURF) 

P-9(NSURF) 


Data for 
aerodynamic 
surface 
no. NSURF 
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P(P)-10(1) 


P(P)-10(NSURF) 
"END" card 


, downwash vector data for 
all surfaces if "LIBR" data 
is not supplied 


Input Format Description for 
Structural and Mode Shape Data Library 

In order to have an input library which contains the struc- 
tural and mode data for all the following aerodynamic problems , 
the first card, L-1, must have the four characters "LIBR" in 
cols. 1-4 on the first card. Anything else may follow on that 
card, for example: 


Cols . 


1-4 

5-80 

LIBR 

ARY FOR F-111 MODES, SWEEP=26 DEG 


It is recommended that the cards be sequenced, however, it is 
not mandatory. 

The second card set, L-2, is two cards in a 6110 format 
which contain the following: 


Cols, (right adjust all data) 6ll0 format 


1-10 

11-20 

21-30 

31-40 

41-50 

51-60 


NSUR 

NMODES 

N1=0 

N2=0 

LI 


L3 

L4 





where 

NF = total number of structural points at which 

mode deflections are given, NF!^ 200, 

NSUR = total number of structural surface spline fits 
to be used for mode interpolation, NSUR<10, 
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NMODES = total number of mode shapes supplied, 

NMODES<20. 

N1=N2=0 = unused constants which must be zero. 

Ll = 0, a comment card of any format will be expected 

at the beginning of each mode shape deck. 

^ 0, no comment card is expected. 

L2 = 0, no mass data will be input for computing 

the generalized masses. 

0, mass data will be input (see Card L-7). 

L3 = 0, no structural points will be modified to 

have zero deflection in each mode. 

0, "L3" total structural points will be modified 
(see Card L-6), L3<NF. 

L4 - 0 , pitch and roll modes will not be added. 

0, (a) the pitch mode will be added as mode 
"NMODES" + 1 with deflections (XF(I), 

1=1, NF). (see Card L-4) 

(b) the roll mode will be added as mode 
"NMODES" + 2 with deflections (YF(I), 

1=1, NF). (see Card L-4) 

The third card set, L-3, is also two cards in a 6F10.0 
format which contains the following: 

Cols. (use decimal or right adjust all data) 6F10.0 format 

Card L-3 


1-10 

11-20 

21-30 

31-40 

41-50 

51-60 

XMODE 

XMASS 

DH 

DWl 

DW2 

BREF 

RHO 
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where 


XMODE 


XMASS 

DH 


EWl 


DW2 


BREF 


RHO 


= uniform mode deflection multiplier. All mode 
deflections are multiplied by this quantity 
which is normally 

XMODE=1.0 

= uniform mass multiplier, normally 

XMASS =1.0 

= uniform multiplier for all deflections 
calculated with the spline interpolation. 

(These are usually the deflections used to 
calculate generalized forces and masses.) 
Normally, DH=1.0 

= uniform multiplier for all slopes calculated 
with the spline interpolation. (These are ” a " 
values in steady flow or the real part of the 
downwash in unsteady flow.) Normally, 

DW1=1.0 

= uniform multiplier for all deflections 
calculated with the spline interpolation 
as the imaginary term in the unsteady 
downwash. Normally, 

DW2=1.0 

= reference length in same dimensions as other 
geometric data. This value must be the same 
as that used in the problem data. 

= density, slugs /ft^. If RHO=0, the default 
value at sea level is used, 

RHO=0. 0023769 slugs/ft^ 

NOTE: for geometric data input in units 

other than inches, a conversion factor 
must be input through either RHO or XMASS, 
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The next card set, L-4, is composed of as many cards as are 
necessary to contain the NF structural point locations in the 
following form; 


read FORMAT 8A10 format 


(input according to FORMAT) 


i XF(1) 

YF(1) i XF(2) I YF(2) [ XF(3) j YF(3) j 

• • • 

XF(NF-l) 

YF(NF-l) XF(NF) 

YF(NF) 



where 


XF(I),YF(I) = X,Y coordinates of the Ith structural point. 

The XF(I) are used to compute the pitch mode for L4?^0 
and YF(I) are used for the roll mode. 

The next card set, L-5, is composed of NSUR subsets, 
each of which contains the following descriptive information 
for each structural surface; 

Cols. (right adjust all data) 6110 and 6F10.0 formats 


Card 

L-5N, 

N=l, 

NSUR 


1-10 

11-20 

21-30 

31-40 41-50 51-60 

NFS (N) 

IF(1,N) 

IF(2,N) 

• • • 

• • • 

IF(NFS(N),N) 

XIS(1,N) XIS(2,N) 

YCS(1,N) X0S(1,N) X0S(2,N) YCS(2,N) 


where 


NFS (N) 


IF(I,N) 

XIS(1,N) 


total number of structural points out of the 
NF set that correspond to the Nth structural 
surface, NFS(N)<100. 

the structural point number assigned to the 
Ith point in the Nth surface. 

X coordinate of the inboard leading edge 
of the Nth structural surface area. 


XIS(2,N) 


X coordinate of the inboard trailing edge. 
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XOS(I,N) = X coordinate of the leading and trailing edges 
at the outboard tip, 

YCS(I,N) = y coordinates of the inboard and outboard 
stations , 

Referring to figure Hi, the spline interpolation scheme 
can be used to fit discontinuous structures. The structural 
points assigned to each surface may be used in one or more 
surfaces. Also, it is not necessary to use all of the NF points. 
The structural surface boundary coordinates XIS, XOS and YCS are 
used to bound the structural point sets. Normally, they are at 
the comers of a square as shovm in the figure. Since the comers 
are transformed to the comers of a square, the boundary coordi- 
nates can also be used to skew or distort the physical coordinate 
system as shown in figure H2, The transformed plane is the plane 
in which the spline fit is made, hence, the decision to distort 
the physical coordinate system should be carefully considered 
since the result may be to introduce errors into the interpolation. 
As an example, if the comer points of a delta wing are used, then 
a simple pitch mode will appear in the 'transformed plane as a 
twisted mode with the pitch angle slope at the root and no pitch 
at the wing tip. Spline interpolation will not usually yield a 
unifom slope distribution in the physical coordinate system for 
such cases. 

If option L3?^0 in card set Lr2, the following data is now 
input. Card set L-6 contains the data for modifying several 
stmctural points such that they have zero deflections for all 
modes with exception of the pitch and roll modes calculated for 
option L4 not zero. The data is input as follows: 

Cols. (right adjust all data) 6110 format 



where 

NF2 number of structural point deflections 

to be given in the mass and mode data 
following this card set, (NOTE: for 

L3=0, NF2=NF) 
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Figure HI.- Structural Surface Subsets 



IF2(I) = structural point which will have zero 
deflection for all modes input. 

This deck permits setting zero deflections at points in addition 
to those input with the mode shapes. As a result, NF2 will 
usually be smaller than NF and the (XF(I) , YF(I)) data given in 
card set L-4 will include the locations of the additional L3 
points as well as the NF2 points. If desired, NF2 can be the 
same as NF and the modal deflections will be zeroed at the IF2(I) 
points in place of the values read as input. 

If option L2?^0, the mass data is read next in card set L-7: 


Card L-7 


where 



AMASS (I) - mass, m^, assigned to the 1^^^ structural point. 

The generalized masses, Mrsj computed as 

NF2 

M„_ = XMASS*XMODE*XMODE*232.0 * V' m.i,. K. 
32.2*RHO*(BREF«3) ^ 


where h^j-, h^g = modal deflections at the i^^ point for 

modes r and s, respectively, 

and XMASS and XMODE are input in card set L-3. 

The deflections are those to be input in the following set. 

Card set L-8 completes the structural library by supplying 
the mode shapes. The form is as follows: 
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(input according to FORMAT on CARD L-7) 



where 


Card L-8 


H(I,J) = modal deflection at the I^h point for the mode. 


Input Format Description for 
Aerod 5 mamic Option and Geometric Data 
Problem Decks 

The problem decks contain all aerod 3 mamic data and the 
structural surfaces used by each aerod 5 mamic surface if the 
library is used. If no library is supplied, the downwash is 
input in the problem decks. In order to distinguish cards that 
are input only when a library is used, these will be denoted as 
P(L)-N where N is the problem card set nximber. For cards used 
only when no library is used, the designation will be P(P)-N. 
Regardless of how many problem decks are used in a single 
computer run, the format which follows will always be used. 

If a library is used, all problems following must be set up 
for use with a library. 

The first card set, P-*!, is a two card title set which is 
printed at the beginning of each problem in the output. The 
format is 

Cols. 

Card P-1 
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Like the library data, it is recommended that the problem decks 
also be sequenced. 

The next card set, P-2, contains the major aerod 5 mamic 
option data as follows : 


Cols , 

(right adjust all data) 

6110 format 

1-10 

11-20 

21-30 31-40 

41-50 

51-60 

lOPl 

I0P2 

I0P3 IOP4 

IQT 

lOPLU 


Card P-2 


where, referring to figure H3, 

IOPl<0, the aerodynamic matrices will be computed 

for all frequencies, 

<0, the matrices will also be written on 

tape unit ” |lOPl| 

>0, the matrices will be read from tape 

unit "lOPl”. 


I0P2 0, aerod 3 mamic results will be obtained in 

the form of pressure distributions, 

^La » ^p» ^cp spanwise distributions 
of lift and chordwise center of pressure, 

=0, generalized forces will be computed for 

"NMODES" total pressure distributions and 
"NMDDES" total integration modes, 

IOP3=IOP4=0, dummy options not used, 

IQT =0, generalized forces will not be written 

on output tape, 

>0, generalized forces will be written on 

tape "IQT" if IOP2=0. 

IOPLU*0, aerodynamic matrices written or read from 

tape "|lOPl|" are in the inverted form, 

?^0, aerod 3 mamic matrices on tape "|lOPl|" 

are in the uninverted form. 
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The next card set, P-3, contains integer data which are 
common to all aerod 5 mamic surfaces. The format is 


3 


where 

NSURF = number of aerod 5 niainic surfaces 

NSURF< 10. 

LS = symmetry option. 

= 0, s)nranetric flow. 

i‘ 0, antis 3 numetric flow, 

NALP = number of alpha or downwash vectors to be input. 
If "LIBR” data is input, NALP=NMDDES< 20. 

NK = number of reduced frequencies. NK 50. 

NW = total number of aerod 3 mamic control points 

and shock load functions for all surfaces 
(see card P-5) 

NW<100, steady flow: NW<70, unsteady flow. 

IDUMP ^ 0 , various intermediate printout is provided 

such as the downwash and integration point 
locations, the uninverted aerodynamic matrices, 
and the pressure series coefficients. 

IDNWSH = 0, downwash vectors are read if "LIBR" data 

is not used. 

^ 0, angle of attack in all downwash vectors 

is set uniformly to ALPO(IS) (given in set P-6N) 
for the ISth surface if "LIBR" data is not used. 


Cols. (right adjust all data) 6110 format 


1 

O 

11-20 

21-30 

31-40 

41-50 

51-60 


NSURF 

LS 

NALP 

NK 

m 

IDUMP 

Card P 

IDNWSH 

JSUROP 
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Figure H3. - Flow Diagram For the Aerodynamic Options 






















JSUROP = 0, regular pressure functions are used 

on all surfaces. 

^ 0, supersonic weighting ftinction is used for 

totally supersonic flow problems. This cannot 
be used for transonic or subsonic problems, 
however, the Mach number distribution may be 
non-uniform as long as no subsonic regions 
are present. 

The next card set contains the real data which are common 
to all aerod}niamic surfaces. The format is 

Cols. (right adjust all data or use decimal) 6F10.0 format 


1-10 

11-20 

21-30 

31-40 41-50 51-60 

XMACH 

BREF 

RK(1) 

• • • 

• • • 

RK(NK) 



where 

XMACH = , free stream Miach number. 

BREF = reference length in same dimensions 

as all geometric data. 

RK(I) = kca = , reduced frequency. 

Uoo 

which completes the data common to all aerod 3 mamic surfaces. 

The following sets are repeated "NSURF” times and they 
begin with card set P-5N given as follows: 


Cols. (right adjust all data) 6110 format 


1-10 

11-20 

21-30 

31-40 

41-50 

51-60 

NC(IS) 

NS (IS) 

NJ(IS) 

NS I (IS) 

ITRANS(IS) 

NLE(IS) 

NTE(IS) 

ICHORD(IS) 

IXI(IS) 

LSPAN(IS) 

LSYM(IS) 

ISTYPE(IS) 

KSURF(IS,1) ... 

• • • 

KSURF(IS,NSURF) 
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where, for the ISth surface, 

NC(IS) * number of chordwise pressure functions 
or downwash points, NC(IS)^10 

NS (IS) = number of spanwise pressure functions 

or downwash chords, and shock load functions 
if KSURF(IS,IS)<0, NS (IS) <15 

I — 

■ ** For steady flow, 

I 

NTR NSURF 

NS(JS) + V NC(IS)*NS(IS) = NW<100 

JS=1 IS=1 

For unsteady flow, 

im NSURF 

23 NS(JS) + 53 NC(IS)*NS(IS) =NW<70 

JS-1 IS=1 

JS » surface for which KSURF(IS ,IS) < 0 (transonic solution) 
NTR “ total number of JS surfaces 


NJ(IS) * number of chordwise integration points: 
NJ(IS) = 0, standard points are used 
NJ(IS) > 0, the input number is used 
as described below for the option IXI(IS). 
NJ(IS)< 15 

NJ(IS) = 0 is recommended for all cases 

NSI(IS)“ number of integration chords on both the 
ISth surface and its image. 

NSI(IS) = 0, standard chords are used 
NSI(IS)>0, the input number is used 
as described in table Hi below for the 
option ISTYPE(IS). 

NSI(IS)< 31 
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ITRANS(IS) 


NLE(IS) 

NTE(IS) 

ICHORD(IS) 

IXI(IS) 


LSPAN(IS) 


LSYM(IS) 


ISTYPE(IS) 


= 0, the Mach nxamber distribution over the 
ISth surface is uniform and eqtial to XMACH 
(card P-4) . 

> 0, the chordwise Mach number distribution 
is read from card P-9N and is used at all 
span stations. 

< 0, the chordwise and spanwise Mach number 
distribution is read from card P-9N. 

* nximber of x,y coordinate pairs to be read 
for defining the leading edge, 

NLE(IS)<10 

= number of points for defining the trailing 
edge. 

NTE(IS)< 10 

= 1,2, 3, 4, the type of chordwise pressure 

distribution to be used for the ISth surface 
as shown in figure H4 . 

= 1,2, 3, 4, the type of chordwise downwash 

point distribution used and since the chord- 
wise integration points are interdigitated, 
this option also determines their distri- 
bution as well. Table H2 summarizes the 
dis tributions . 


IXI(IS)=1 is recommended 
for almost all cases 


= 1,2, 3, 4, the type of spanwise loading 
to be used for the ISth surface as shown 
in figure H4. 

= 0, S 3 nnmetric flow. 

+ 0, antisymmetric flow. 

= 1,2,3, the type of surface and image 
to be usedo Type 1 has continuous 
functions across the centerline whereas 
types 2 and 3 are independent. The 
type 2 spanwise loading can be used with 
all three surface types. 
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TABLE H2. CHORDWISE DOWNWASH AND INTEGRATION POINTS (m=NC(IS)) 
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For supersonic flow the are not constant and are a function of the downwash 
point location. For interference, the ^ are constant as determined by IXI 
and NJ default is NJ*2*NC(IS), If NJ is input, then none of the restrictions 
given in the table for NJ=J apply. 
















TYPE 1 TYPE 2 TYPE 3 


SURFACE TYPES 
ISTYPE (IS) 


tin) h(o 



SPANWISE LOADINGS CHORDWISE LOADINGS 

LSPAN (IS) ICHORD (IS) 


Figure H4, - Surface and Loading Types in the 
Aerodynamic Program 
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NS I 
















KSURF(IS,KS) = interference calculation option as described 

below and illustrated in figure H5. 

For IS = KS: 

KSURF(IS,IS) = 0, surface IS is not a transonic surface. 

KSURF(IS,IS) > 0, surface IS is the upstream svirface in a 

transonic pair with a shock along its 
trailing edge. The other surface in the 
pair is the value of KSURF(IS ,IS) . 

(see figure H5) 

KSURF(IS ,IS) < 0, surface IS is the downstream surface 

in a transonic pair with a shock at its 
leading edge. The other surface is the 
absolute value of KSURF(IS ,IS) . 

(see figure H5) 

For IS KS, 

KSURF(IS,KS) =* 0, interference of surface KS on IS will be 

calculated but spanwise integration correc- 
tion terms will not be calculated. This 
option is usually taken if KS is downstream 
of or parallel to IS. 

< 0, no interference is calculated for 
KS on IS. 

> 0, interference is calculated with the 
spanwise integral correction terms. 

The next card set in the NSURF decks is P-6N which contains 
the following real data: 

Cols. (right adjust all data or use decimal) 6F10.0 format 


1-10 

11-20 

21-30 

31-40 

41-50 

51-60 

ALPO(IS) 

XV(IS) 

YV(IS) 

ZV(IS) 

THETA (IS) 
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Figure H5. KSURF(IS,KS) For Transonic Flow Over 
Interfering Surfaces 


L23 




where 


ALPO(IS) 


XV(IS) 

YV(IS) 

ZV(IS) 

A 


- uniform angle of attack, in radians, 
which is always added to all downwash 
vectors for all control points that fall 
on the ISth surface. 

^ x,y,z coordinates of the leading edge 
“ inboard tip of the surface (see figure H6) , 


THETA (IS) = inclination of the surface relative to 

the y-axis (see figure H6) , 


The planform geometry data is given next in card set P-7N 
as follows : 


Cols. (right adjust all data or use decimal) 6F10.0 format 


1-10 

11-20 21-30 
1 

31-40 

41-50 51-60 

XLE(l) 

YLE(l) 

XLE (2) 

YLE(2). 

• • • 

• • • 

XLE(NLE) 

YLE(NLE) 


XTE(l) 

YTE(l) 

XTE (2) 

YTE (2) 

• • • 

• • • 

XTE(NTE) 

YTE(NTE) 



where 


XLE(I) , YLE(I) = x,y coordinates of the Ith leading edge 

break point. It is not necessary for 
XLE(1)=XV(IS) and YLE(1)=W(IS) . 

XTE(I) , YTE(I) = x,y coordinates of the Ith trailing 

edge break point » 

The dimensions of XLE, YLE, XTE, YTE are the same as all other 
input geometric data. 


If the structural and mode shape library is used, the 
following card must be supplied: 


Cols. (right adjust) IlO format 


1-10 


NST(IS) 



Card P(L)-8N 
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XV(IS) 
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Figure H6. - Aerodynamic Geometry Input Data For Arbitrary Surfaces 



where 


NST(IS) = the structural surface given in the library 
from which deflections and downwash are 
calculated for the mode shapes. 

If library data is not supplied then this card must be omitted. 

If ITRANS(IS)?^0, then the Mach number distribution data is 
supplied in card set P-9N, For ITRANS(IS) > 0, the chordwise 
distribution is assumed constant along the span and the data is 
given as 


Cols, (right adjust or use decimal) 6F10.0 format 


1-10 

11-20 

21-30 

31-40 

41-50 

51-60 

GMACH(1,1) 

GMACH(2,1) 

• • o 




• • • 

GMACH(Mc,l) 





Card P-9N 
(ITRANS > 0) 

where 


GMACH(I,1) = input Mach nvimber at constant percent 
chord line /I-l s 
^Mc-1'' 

Me = NC(IS) + 1 

Note that the data is input at constant chord fraction increments 
from the leading to the trailing edge. For ITRANS (IS) < 0, the 
distribution is input at (NS(IS)+1) span stations as follows: 


Cols. (right adjust or use decimal) 6F10,0 format 


1-10 

11-20 

21-30 

31-40 

41-50 

51-60 

GMACH(1,1) 

GMACH(2,1) 

• • • 

GMACH(Mc,l) 


• • • 

GMACH(1,MS) 

GMACH(2,MS) 

• • • 

• • • 

GMACH(Mc,MS) 



Card P-9N 
(ITRANS < 0) 
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where 


GMACH(I,J) = input number at the (~ - j, ) 


and spanwise locations. 


chordwise 


MS = NS (IS) + 1 

Again, the spanwise distribution is uniformly spaced starting 
from the inboardmost station to and including the tip at / 1 n 
increments. ''MS-1'^ 

This completes the input aerod 3 mamic geometry and option 
data and the problem data if the structural library is used. 

If the library is not used then the downwash must be input after 
all aerodynamic data for all surfaces, if the IDNWSH option in 
card P-3 is not used. For steady flow, the data are input as 
follows : 


Cols. (right adjust or use decimal) 6F10.0 format 


1-10 

11-20 

21-30 

31-40 

41-50 

51-60 

ALP(1,1,IS) 

ALP(2,1,IS) 

• • • 

• • • 

ALP(NI,1,IS) 


... 

^LP(1,NALP,IS) 

ALP (2, NALP, IS) 

• • • 

• 0 • 

ALP (NI, NALP, IS) 


where 


Card P(P)-10NN 
k=0 


ALP(I,K,IS) = slope (radians) at the Ith downwash point 
in the Kth downwash vector on the ISth 
surface, from LE to TE , outbd to inbd . 

NALP = total niamber of downwash vectors as input 

on card P-3, NALP <20. 

NI = NC(IS)*NS (IS) , total number of downwash 

points on the ISth surface. 


127 




For unsteady flow, the only change is the addition of deflections: 


Cols. (right adjust or use decimal) 6F10.0 format 


1-10 

11-20 

21-30 

31-40 

41-50 

51-60 

ALP(1,1,IS) 

ALP (2,1, IS) 

• • • 

• • • 

ALP(NI,1,IS) 


H(1,1,IS) 

H(2,1,IS) 

• • • 

• • 

H(NI,1,IS) 


1 

• • • 

1 

ALP(1,NALP,IS) 

ALP(2,NALP,IS) 

• • • 

• • • 

ALP(NI,NALP,IS) 


H(1,NALP,IS) 

H(2,NALP,IS) 

• • • 

• • • 

H(NI,NALP,IS) 



Card P(P)-10NN 

where k?^0 


H(I,K,IS) * deflection (same dimensions as the 
geometric data) at the Ith downwash 
point in the Kth downwash vector on 
the ISth surface. 

In both cases, the downwash data are given only after all of the 
aerod 3 mamic geometry and option data are completed for all 
surfaces. Then, as the format shows, all modal data is input 
for each surface together. 

The last card in a problem deck signals the program to 
continue to the next problem or terminate the job. This is 
the "END" card which terminates the job. 


Card 

"END" 


Cols , 


1-4 

5-60 

ENDb 
1 , 

(JOB WILL BE TERMINATED) 


where b is a blank space. For continuation, the characters 
in Cols. 1-4 may be anything but ENDb, for example 


Cols , 


1-4 

5-60 

PROB 

LEM 1 END (CONT. ON TO NEXT PROBLEM) 


Card 

"END" 


which will continue the job. 
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